Paso 5 (1)- Trayectorias de hospitalización y mortalidad con foco en condiciones vinculadas a trastornos de salud mental y consumo de sustancias posterior a un primer ingreso por alguno de estos trastornos, en usuarios/as jóvenes y adultos emergentes de población general y pertenecientes a pueblos originarios, 2018-2021, Chile
Representar las mejores opciones de agrupamiento, junto con su relación con otras variables. En esta etapa, analizamos la asociación con covariables de la solución con mejores índices de calidad para la resolución trimestral.
Autor/a
Andrés González Santa Cruz
Fecha de publicación
25 de ene, 2025
Configurar
Código
# remover objetos y memoria utilizadarm(list=ls());gc()
used (Mb) gc trigger (Mb) max used (Mb)
Ncells 600238 32.1 1291560 69 707674 37.8
Vcells 1163664 8.9 8388608 64 1876475 14.4
#elegir repositorioif(Sys.info()["sysname"]=="Windows"){options(repos =c(CRAN ="https://cran.dcc.uchile.cl/"))}options(install.packages.check.source ="yes") # Chequea la fuente de los paquetes#borrar caché#system("fc-cache -f -v")if(!require(pacman)){install.packages("pacman");require(pacman)}pacman::p_unlock(lib.loc =.libPaths()) #para no tener problemas reinstalando paquetesif(Sys.info()["sysname"]=="Windows"){if (getRversion() !="4.4.0") { stop("Requiere versión de R 4.4.0. Actual: ", getRversion()) }}if(!require(job)){install.packages("job");require(job)}if(!require(kableExtra)){install.packages("kableExtra");require(kableExtra)}if(!require(tidyverse)){install.packages("tidyverse");require(tidyverse)}if(!require(cluster)){install.packages("cluster"); require(cluster)}if(!require(WeightedCluster)){install.packages("WeightedCluster"); require(WeightedCluster)}if(!require(devtools)){install.packages("devtools"); require(devtools)}if(!require(TraMineR)){install.packages("TraMineR"); require(TraMineR)}if(!require(TraMineRextras)){install.packages("TraMineRextras"); require(TraMineRextras)}if(!require(NbClust)){install.packages("NbClust"); require(NbClust)}if(!require(haven)){install.packages("haven"); require(haven)}if(!require(ggseqplot)){install.packages("ggseqplot"); require(ggseqplot)}if(!require(grid)){install.packages("grid"); require(grid)}if(!require(gridExtra)){install.packages("gridExtra"); require(gridExtra)}if(!require(Tmisc)){install.packages("Tmisc"); require(Tmisc)}if(!require(factoextra)){install.packages("factoextra"); require(factoextra)}if(!require(stargazer)){install.packages("stargazer"); require(stargazer)}if(!require(gtsummary)){install.packages("gtsummary"); require(gtsummary)}if(!require(lmtest)){install.packages("lmtest"); require(lmtest)}if(!require(emmeans)){install.packages("emmeans"); require(emmeans)}if(!require(fpp2)){install.packages("fpp2"); require(fpp2)}if(!require(purrr)){install.packages("purrr"); require(purrr)}if(!require(forecast)){install.packages("forecast"); require(forecast)}if(!require(magrittr)){install.packages("magrittr"); require(magrittr)}if(!require(foreach)){install.packages("foreach"); require(foreach)}if(!require(doParallel)){install.packages("doParallel"); require(doParallel)}if(!require(progressr)){install.packages("progressr"); require(progressr)}if(!require(chisq.posthoc.test)){devtools::install_github("ebbertd/chisq.posthoc.test")}if(!require(rstatix)){install.packages("rstatix"); require(rstatix)}if(!require(rio)){install.packages("rio"); require(rio)}if(!require(cowplot)){install.packages("cowplot"); require(cowplot)}if(!require(DiagrammeR)){install.packages("DiagrammeR"); require(DiagrammeR)}if(!require(DiagrammeRsvg)){install.packages("DiagrammeRsvg"); require(DiagrammeRsvg)}if(!require(rsvg)){install.packages("rsvg"); require(rsvg)}if(!require(survminer)){install.packages("survminer"); require(survminer)}seq_mean_t_dos_grupos <-function(bd =NULL, group1, group2) {# Agrupar por ambas variables resultados <-by(bd, list(group1, group2), seqmeant)# Obtener todas las combinaciones posibles de los grupos combinaciones <-expand.grid(group1 =unique(group1), group2 =unique(group2), stringsAsFactors =FALSE)# Extraer los resultados y asociarlos con las combinaciones resultados_df <-do.call(rbind, lapply(seq_along(resultados), function(i) { group_name1 <-attr(resultados, "dimnames")[[1]][i] group_name2 <-attr(resultados, "dimnames")[[2]][i]data.frame(factor_inclusivo_1 = group_name1, factor_inclusivo_2 = group_name2, Mean = resultados[[i]]) }))# Unir los resultados con las combinaciones para rellenar los valores faltantes final_df <-merge(combinaciones, resultados_df, by.x =c("group1", "group2"), by.y =c("factor_inclusivo_1", "factor_inclusivo_2"), all.x =TRUE)return(final_df)}multinom_pivot_wider <-function(x) {# check inputs match expectatations# create tibble of results df <- tibble::tibble(outcome_level =unique(x$table_body$groupname_col)) df$tbl <- purrr::map( df$outcome_level,function(lvl) { gtsummary::modify_table_body( x, ~dplyr::filter(.x, .data$groupname_col %in% lvl) %>% dplyr::ungroup() %>% dplyr::select(-.data$groupname_col) ) } )tbl_merge(df$tbl, tab_spanner =paste0("**", df$outcome_level, "**"))}best_subset_multinom <-function(y, x.vars, data) {# y Nombre de la variable dependiente (cadena de texto)# x.vars Vector de nombres de predictores (caracter)# data Dataframe con los datos de entrenamiento# Cargar las librerías necesariasrequire(dplyr)require(purrr)require(tidyr)require(nnet)require(MASS)# Generar todas las combinaciones posibles de predictores predictors_list <-lapply(1:length(x.vars), function(i) {combn(x.vars, i, simplify =FALSE) }) %>%unlist(recursive =FALSE)# Inicializar una lista para almacenar los resultados results <-list()# Iterar sobre cada combinación de predictoresfor (i inseq_along(predictors_list)) { predictors <- predictors_list[[i]] formula <-as.formula(paste(y, "~", paste(predictors, collapse ="+")))# Ajustar el modelo multinomial model <-tryCatch( nnet::multinom(formula, data = data, trace =FALSE),error =function(e) NULL )# Si el modelo se ajustó correctamente, almacenar los resultadosif (!is.null(model)) {# Extraer el AIC del modelo aic <-AIC(model)# Almacenar la información en una lista results[[length(results) +1]] <-list(predictors = predictors,model = model,AIC = aic ) } }# Convertir la lista de resultados en un dataframe results_df <- results %>% purrr::map_df(function(res) {data.frame(predictors =paste(res$predictors, collapse ="+"),AIC = res$AIC,stringsAsFactors =FALSE ) })# Ordenar los modelos por AIC de menor a mayor results_df <- results_df %>%arrange(AIC)return(results_df)}best_subset_multinom_interactions <-function(y, x.vars, data) {# y Nombre de la variable dependiente (cadena de texto)# x.vars Vector de nombres de predictores (caracter)# data Dataframe con los datos de entrenamiento# Cargar las librerías necesariasrequire(dplyr)require(purrr)require(tidyr)require(nnet)require(MASS)# Generar todas las combinaciones posibles de predictores (efectos principales) main_effects_list <-lapply(1:length(x.vars), function(i) {combn(x.vars, i, simplify =FALSE) }) %>%unlist(recursive =FALSE)# Inicializar una lista para almacenar los resultados results <-list()# Iterar sobre cada combinación de efectos principalesfor (main_effects in main_effects_list) {# Generar términos de interacción de hasta 3 variables interaction_terms <-list()# Para interacciones de 2 variablesif (length(main_effects) >=2) { interaction_terms_2way <-combn(main_effects, 2, function(x) paste(x, collapse =":")) interaction_terms <-c(interaction_terms, interaction_terms_2way) }# Para interacciones de 3 variablesif (length(main_effects) >=3) { interaction_terms_3way <-combn(main_effects, 3, function(x) paste(x, collapse =":")) interaction_terms <-c(interaction_terms, interaction_terms_3way) }# Combinar efectos principales e interacciones all_terms <-c(main_effects, interaction_terms)# Generar todas las combinaciones posibles de términos (incluyendo interacciones)# Solo se incluyen interacciones si sus efectos principales están presentes term_combinations <-list()# Obtener todos los subconjuntos de efectos principales main_effects_subsets <-lapply(1:length(main_effects), function(i) {combn(main_effects, i, simplify =FALSE) }) %>%unlist(recursive =FALSE)# Para cada subconjunto de efectos principalesfor (me in main_effects_subsets) {# Iniciar con los efectos principales terms <- me# Incluir interacciones solo si todos sus efectos principales están incluidos possible_interactions <- interaction_terms[sapply(interaction_terms, function(x) { vars_in_interaction <-unlist(strsplit(x, ":"))all(vars_in_interaction %in% me) }) ]# Generar todas las combinaciones de interacciones para incluir interaction_subsets <-list(NULL)if (length(possible_interactions) >0) { interaction_subsets <-lapply(1:length(possible_interactions), function(i) {combn(possible_interactions, i, simplify =FALSE) }) %>%unlist(recursive =FALSE) }# Para cada combinación de interacciones, crear el conjunto completo de términosfor (ints in interaction_subsets) {if (is.null(ints)) { full_terms <- terms } else { full_terms <-c(terms, ints) }# Añadir a la lista de combinaciones de términos term_combinations <-append(term_combinations, list(full_terms)) } }# Ajustar modelos para cada combinación de términosfor (terms in term_combinations) { formula <-as.formula(paste(y, "~", paste(terms, collapse ="+")))# Ajustar el modelo multinomial model <-tryCatch( nnet::multinom(formula, data = data, trace =FALSE),error =function(e) NULL,warning =function(w) NULL )# Si el modelo se ajustó correctamente, almacenar los resultadosif (!is.null(model)) {# Extraer el BIC del modelo bic <-BIC(model)# Almacenar la información en la lista de resultados results[[length(results) +1]] <-list(predictors =paste(terms, collapse =" + "),model = model,BIC = bic ) } } }# Convertir la lista de resultados en un dataframe results_df <- results %>% purrr::map_df(function(res) {data.frame(predictors = res$predictors,BIC = res$BIC,stringsAsFactors =FALSE ) })# Ordenar los modelos por BIC de menor a mayor results_df <- results_df %>%arrange(BIC)return(results_df)}best_subset_multinom_interactions_parallel <-function(y, x.vars, data) {# y Nombre de la variable dependiente (cadena de texto)# x.vars Vector de nombres de predictores (caracter)# data Dataframe con los datos de entrenamiento# Cargar las librerías necesarias dentro de la funciónrequire(dplyr)require(purrr)require(tidyr)require(nnet)require(MASS)require(foreach)require(doParallel)require(progressr)# Iniciar los gestores de progresohandlers(global =TRUE)handlers("txt")# Generar todas las combinaciones posibles de predictores (efectos principales) main_effects_list <-lapply(1:length(x.vars), function(i) {combn(x.vars, i, simplify =FALSE) }) %>%unlist(recursive =FALSE)# Inicializar una lista para almacenar las fórmulas de los modelos formulas_list <-list()# Generar todas las fórmulas posibles con interacciones hasta de 3 variablesfor (main_effects in main_effects_list) {# Generar términos de interacción de hasta 3 variables interaction_terms <-character(0) # Aseguramos que es un vector de caracteres# Para interacciones de 2 variablesif (length(main_effects) >=2) { interaction_terms_2way <-combn(main_effects, 2, function(x) paste(x, collapse =":"), simplify =TRUE) interaction_terms <-c(interaction_terms, interaction_terms_2way) }# Para interacciones de 3 variablesif (length(main_effects) >=3) { interaction_terms_3way <-combn(main_effects, 3, function(x) paste(x, collapse =":"), simplify =TRUE) interaction_terms <-c(interaction_terms, interaction_terms_3way) }# Generar todas las combinaciones posibles de efectos principales main_effects_subsets <-lapply(1:length(main_effects), function(i) {combn(main_effects, i, simplify =FALSE) }) %>%unlist(recursive =FALSE)# Para cada subconjunto de efectos principalesfor (me in main_effects_subsets) {# Iniciar con los efectos principales terms <- me# Identificar interacciones cuyos efectos principales están en 'me'if (length(interaction_terms) >0) { possible_interactions <- interaction_terms[vapply(interaction_terms, function(x) { vars_in_interaction <-unlist(strsplit(x, ":"))all(vars_in_interaction %in% me) }, FUN.VALUE =logical(1)) ] } else { possible_interactions <-character(0) }# Generar todas las combinaciones posibles de estas interacciones interaction_subsets <-list(character(0)) # Incluir el caso sin interaccionesif (length(possible_interactions) >0) { interaction_combinations <-lapply(1:length(possible_interactions), function(i) {combn(possible_interactions, i, simplify =FALSE) }) %>%unlist(recursive =FALSE) interaction_subsets <-c(interaction_subsets, interaction_combinations) }# Para cada combinación de interaccionesfor (ints in interaction_subsets) { full_terms <-c(terms, ints)# Crear la fórmula del modelo y almacenarla formula_str <-paste(y, "~", paste(full_terms, collapse ="+")) formulas_list <-append(formulas_list, list(formula_str)) } } }# Eliminar posibles duplicados de fórmulas formulas_list <-unique(formulas_list)# Total de modelos a ajustar total_models <-length(formulas_list)# Iniciar el progreso p <-progressor(steps = total_models)# Ajustar los modelos en paralelo usando foreach results_list <-foreach(i =1:total_models, .packages =c("nnet", "MASS"), .combine ='rbind') %dopar% { formula_str <- formulas_list[[i]] formula <-as.formula(formula_str)# Ajustar el modelo model <-tryCatch( nnet::multinom(formula, data = data, trace =FALSE),error =function(e) NULL,warning =function(w) NULL )# Actualizar el progresop(sprintf("Ajustando modelo %d de %d", i, total_models))# Si el modelo se ajustó correctamente, almacenar los resultadosif (!is.null(model)) { bic <-BIC(model)data.frame(predictors = formula_str,BIC = bic,stringsAsFactors =FALSE ) } else {NULL } }# Convertir los resultados a dataframe y ordenar por BIC results_df <-as.data.frame(results_list) results_df <- results_df %>%arrange(BIC)return(results_df)}num_cores <- parallel::detectCores() -1cl <-makeCluster(num_cores)registerDoParallel(cl)#pacman job kableExtra tidyverse cluster WeightedCluster devtools TraMineR TraMineRextras NbClust haven ggseqplot gridExtra Tmisc factoextra reticulate withr rmarkdown quartooptions(knitr.kable.NA ='')#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:##:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#knitr::knit_hooks$set(time_it =local({ now <-NULLfunction(before, options) {if (before) {# record the current time before each chunk now <<-Sys.time() } else {# calculate the time difference after a chunk res <-ifelse(difftime(Sys.time(), now)>(60^2),difftime(Sys.time(), now)/(60^2),difftime(Sys.time(), now)/(60^1))# return a character string to show the time x<-ifelse(difftime(Sys.time(), now)>(60^2),paste("Tiempo que demora esta sección:", round(res,1), "horas"),paste("Tiempo que demora esta sección:", round(res,1), "minutos"))paste('<div class="message">', gsub('##', '\n', x),'</div>', sep ='\n') } }}))knitr::opts_chunk$set(time_it =TRUE)#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:format_cells <-function(df, rows ,cols, value =c("italics", "bold", "strikethrough")){# select the correct markup# one * for italics, two ** for bold map <-setNames(c("*", "**", "~~"), c("italics", "bold", "strikethrough")) markup <- map[value] for (r in rows){for(c in cols){# Make sure values are not factors df[[c]] <-as.character( df[[c]])# Update formatting df[r, c] <-ifelse(nchar(df[r, c])==0,"",paste0(markup, gsub(" ", "", df[r, c]), markup)) } }return(df)}#To produce line breaks in messages and warningsknitr::knit_hooks$set(error =function(x, options) {paste('\n\n<div class="alert alert-danger">',gsub('##', '\n', gsub('^##\ Error', '**Error**', x)),'</div>', sep ='\n') },warning =function(x, options) {paste('\n\n<div class="alert alert-warning">',gsub('##', '\n', gsub('^##\ Warning:', '**Warning**', x)),'</div>', sep ='\n') },message =function(x, options) {paste('<div class="message">',gsub('##', '\n', x),'</div>', sep ='\n') })#_#_#_#_#_#_#_#_#_#_#_#_#_invisible("Function to format CreateTableOne into a database")as.data.frame.TableOne <-function(x, ...) {capture.output(print(x,showAllLevels =TRUE, varLabels = T,...) -> x) y <-as.data.frame(x) y$characteristic <- dplyr::na_if(rownames(x), "") y <- y %>%fill(characteristic, .direction ="down") %>% dplyr::select(characteristic, everything())rownames(y) <-NULL y}#_#_#_#_#_#_#_#_#_#_#_#_#_# Austin, P. C. (2009). The Relative Ability of Different Propensity # Score Methods to Balance Measured Covariates Between # Treated and Untreated Subjects in Observational Studies. Medical # Decision Making. https://doi.org/10.1177/0272989X09341755smd_bin <-function(x,y){ z <- x*(1-x) t <- y*(1-y) k <-sum(z,t) l <- k/2return((x-y)/sqrt(l))}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:if(.Platform$OS.type =="windows") withAutoprint({memory.size()memory.size(TRUE)memory.limit()})
> memory.size()
[1] Inf
> memory.size(TRUE)
[1] Inf
> memory.limit()
[1] Inf
Código
#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:func_tab_range_clus<-function(range_clus){rbind.data.frame(lapply(list(as.vector(rev(sort(table(range_clus$clustering$cluster2)))),as.vector(rev(sort(table(range_clus$clustering$cluster3)))),as.vector(rev(sort(table(range_clus$clustering$cluster4)))),as.vector(rev(sort(table(range_clus$clustering$cluster5)))),as.vector(rev(sort(table(range_clus$clustering$cluster6)))),as.vector(rev(sort(table(range_clus$clustering$cluster7)))),as.vector(rev(sort(table(range_clus$clustering$cluster8)))),as.vector(rev(sort(table(range_clus$clustering$cluster9)))),as.vector(rev(sort(table(range_clus$clustering$cluster10)))),as.vector(rev(sort(table(range_clus$clustering$cluster11)))),as.vector(rev(sort(table(range_clus$clustering$cluster12)))),as.vector(rev(sort(table(range_clus$clustering$cluster13)))),as.vector(rev(sort(table(range_clus$clustering$cluster14)))),as.vector(rev(sort(table(range_clus$clustering$cluster15)))) ),function(x) { length_out <-max(sapply(list(as.vector(rev(sort(table(range_clus$clustering$cluster2)))),as.vector(rev(sort(table(range_clus$clustering$cluster3)))),as.vector(rev(sort(table(range_clus$clustering$cluster4)))),as.vector(rev(sort(table(range_clus$clustering$cluster5)))),as.vector(rev(sort(table(range_clus$clustering$cluster6)))),as.vector(rev(sort(table(range_clus$clustering$cluster7)))),as.vector(rev(sort(table(range_clus$clustering$cluster8)))),as.vector(rev(sort(table(range_clus$clustering$cluster9)))),as.vector(rev(sort(table(range_clus$clustering$cluster10)))),as.vector(rev(sort(table(range_clus$clustering$cluster11)))),as.vector(rev(sort(table(range_clus$clustering$cluster12)))),as.vector(rev(sort(table(range_clus$clustering$cluster13)))),as.vector(rev(sort(table(range_clus$clustering$cluster14)))),as.vector(rev(sort(table(range_clus$clustering$cluster15)))) ), length))c(x, rep(NA, length_out -length(x))) } ))%>%t() |>data.frame()%>%`rownames<-`(NULL)}frobenius_norm <-function(matrix1, matrix2) {if (!all(dim(matrix1) ==dim(matrix2))) {stop("Matrices must have the same dimensions") }# Replace NA values with 0 (or any other desired default) matrix1[is.na(matrix1)] <-0 matrix2[is.na(matrix2)] <-0# Calculate the residuals residuals <- matrix1 - matrix2# Frobenius norm frobenius <-sqrt(sum(residuals^2))return(frobenius)}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:confcqi2 <-function(nullstat, quant, n){ alpha <- (1-quant)/2#calpha <- alpha+(alpha-1)/n#print(c(calpha, alpha))#minmax <- quantile(nullstat, c(calpha, 1-calpha)) minmax <-quantile(nullstat, c(alpha, 1-alpha))return(minmax)}normstatcqi2 <-function(bcq, stat, norm=TRUE){ origstat <- bcq$clustrange$stats[, stat] nullstat <- bcq$stats[[stat]]#normstat <- rbind(nullstat, origstat)if(norm){for(i inseq_along(origstat)){ mx <-mean(nullstat[, i]) sdx <-sd(nullstat[, i]) nullstat[ , i] <- (nullstat[, i]-mx)/sdx origstat[i] <- (origstat[i]-mx)/sdx } } alldatamax <-apply(nullstat, 1, max)#as.vector(xx) sumcqi <-list(origstat=origstat, nullstat=nullstat, alldatamax=alldatamax)return(sumcqi)}print.seqnullcqi.powder <-function(x, norm =FALSE, quant =0.95, digits =2, append =FALSE, ...) {cat("Parametric bootstrap cluster analysis validation\n")cat("Sequence analysis null model:", deparse(x$nullmodel), "\n")cat("Number of bootstraps:", x$R, "\n")cat("Clustering method:", ifelse(x$kmedoid, "PAM/K-Medoid", paste0("hclust with ", x$hclust.method)), "\n")cat("Seqdist arguments:", deparse(x$seqdist.args), "\n\n\n") alls <-as.data.frame(x$clustrange$stats) quants <-rep("", ncol(alls))names(quants) <-colnames(alls)for (ss incolnames(alls)) { sumcqi <-normstatcqi2(x, stat = ss, norm = norm) alls[, ss] <-as.character(round(sumcqi$origstat, digits = digits)) borne <-as.character(round(confcqi2(sumcqi$alldatamax, quant, x$R), digits = digits)) quants[ss] <-paste0("[", borne[1], "; ", borne[2], "]") } results_tibble <- tibble::as_tibble(rbind(alls, rep("", length(quants)), quants))# Print a summary to the console for immediate feedbackrownames(results_tibble) <-c(rownames(x$clustrange$stats), "", paste("Null Max-T", quant, "interval")) results_df <-as.data.frame(results_tibble)print(results_tibble, ...)return(list(results_tibble= results_tibble, results_df= results_df ))}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:##:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:# Función para aplicar la prueba de Fisher a todas las combinaciones de filas usando todas las columnasfisher_posthoc_all_cols <-function(contingency_table) {# Obtener combinaciones de filas (pares) row_pairs <-combn(rownames(contingency_table), 2, simplify =FALSE)# Aplicar la prueba de Fisher a cada par de filas usando todas las columnas al mismo tiempo results <-map_dfr(row_pairs, function(pair) {# Crear tabla de 2xN para el par de filas en todas las columnas sub_table <- contingency_table[pair, , drop =FALSE]# Aplicar el test de Fisher test_result <-fisher.test(sub_table, simulate.p.value=T,B=1e4)# Devolver los resultados en un data frametibble(Row1 = pair[1],Row2 = pair[2],p.value = test_result$p.value ) })# Ajustar p-valores usando el método de Holm results <- results %>%mutate(p.adjusted =p.adjust(p.value, method ="holm"))return(results)}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:save_base_plot_as_grob <-function(plot_expr, res=300, width =1600, height=1200) {# Crea un archivo temporal con extensión .png filename <-tempfile(fileext =".png")# Guarda el gráfico en alta resolución en el archivo temporalpng(filename, width = width, height = height, res = res)replayPlot(plot_expr) # Reproduce el gráfico grabadodev.off() # Cierra el dispositivo gráfico# Convierte el archivo PNG en un objeto gráfico (grob) grob <- grid::rasterGrob(png::readPNG(filename), interpolate =TRUE)return(grob) # Devuelve el grob}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:chisq_cramerv<-function(contingency_table){ chisq_test <-chisq.test(contingency_table) cramers_v <-sqrt(chisq_test$statistic / (sum(contingency_table) * (min(dim(contingency_table)) -1)))list(chisq_statistic=sprintf("%1.2f", chisq_test$statistic), chisq_df= chisq_test$parameter, chisq_p_value =ifelse(chisq_test$p.value<.001, "<0.001", sprintf("%1.4f", chisq_test$p.value)), cramers_v =sprintf("%1.2f", cramers_v))}#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#oneway_anova_effect_size <-function(values, group) {# Perform one-way ANOVA anova_result <-aov(values ~ group)# Summarize ANOVA results anova_summary <-summary(anova_result)# Extract sums of squares ss_between <- anova_summary[[1]]$"Sum Sq"[1] ss_total <-sum(anova_summary[[1]]$"Sum Sq")# Calculate eta-squared eta_squared <- ss_between / ss_total# Return ANOVA summary and effect sizelist(anova_summary = anova_summary,eta_squared = eta_squared )}
Resultados
0.a. Historial
A continuación, mostramos un ejemplo para utilizar el algoritmo LCS (Longest common subsequence o Subsecuencia común más larga).
Código
invisible("Soluciones de cluster seleccionadas")# States_Wide.seq_quarter_t_prim_adm_cens, # group=om_dist_quarter_c$clustering$cluster9, # # # tates_Wide.seq_quarter_t_prim_adm_cens, # group=om_dist_quarter_c$clustering$cluster6 # 2do lugar# # States_Wide.seq_quarter_t_prim_adm_cens, # group=pamRange_quarter_om$clustering$cluster6# # # States_Wide.seq_month_t_prim_adm_cens, # group=lcs_dist_month_c$clustering$cluster4 ## # # States_Wide.seq_quarter_t_prim_adm_cens, # group=lcs_dist_quarter_c_rm$clustering$cluster5# # # pamRange_quarter_om$clustering$cluster6# # States_Wide.seq_quarter_t_prim_adm_RM# invisible("Bases sobre las que utilizar")#dt_ing_calendar_quarter_t_desde_primera_adm_dedup##dt_ing_calendar_month_t_desde_primera_adm_dedup## Definir las cadenas ficticiascadena_A <-c("cp", "aus", "aus", "psi", "psi", "cens")cadena_B <-c("cp", "aus", "psi", "psi", "psi", "aus")# Longitudes de las cadenasn <-length(cadena_A)m <-length(cadena_B)# Crear la matriz LCSLCS <-matrix(0, nrow = n +1, ncol = m +1)# Llenar la matriz según el algoritmo LCSfor (i in2:(n +1)) {for (j in2:(m +1)) {if (cadena_A[i -1] == cadena_B[j -1]) { LCS[i, j] <- LCS[i -1, j -1] +1 } else { LCS[i, j] <-max(LCS[i -1, j], LCS[i, j -1]) } }}# Convertir la matriz a un data frame para ggplotLCS_df <-as.data.frame(LCS) %>%mutate(row =0:n) %>%pivot_longer(-row, names_to ="col", values_to ="value") %>%mutate(col =as.integer(gsub("V", "", col)) -1)# Etiquetas para las cadenaslabels_A <-c("", cadena_A)labels_B <-c("", cadena_B)# Graficar con ggplotggplot(LCS_df, aes(x = col, y = row, fill = value)) +geom_tile(color ="white") +geom_text(aes(label = value), color ="black", size =6) +scale_fill_gradient(low ="white", high ="blue") +scale_x_continuous(breaks =0:m, labels = labels_B) +scale_y_reverse(breaks =0:n, labels = labels_A) +labs(title ="Matriz de LCS para cadenas ['cp', 'aus', 'aus', 'psi', 'psi', 'cens'] y\n ['cp', 'aus', 'psi', 'psi', 'psi', 'aus']",fill ="Subsecuencia\ncomún más\nlarga",x ="Cadena B",y ="Cadena A") +theme_minimal()+theme(plot.title =element_text(size =20, face ="bold"),axis.title =element_text(size =16),axis.text =element_text(size =14),legend.title =element_text(size =14),legend.text =element_text(size =12) )
Ejemplo LCS
Tiempo que demora esta sección: 0 minutos
De igual forma, intentamos graficar
Código
import numpy as npimport matplotlib.pyplot as pltfrom sklearn.metrics import pairwise_distancesfrom sklearn.cluster import KMeansimport seaborn as snsfrom matplotlib.patches import Ellipse# Generar datos artificiales bien centradosnp.random.seed(42)cluster_1 = np.random.normal(loc=0, scale=0.5, size=(30, 2))cluster_2 = np.random.normal(loc=5, scale=0.5, size=(30, 2))cluster_3 = np.random.normal(loc=10, scale=0.5, size=(30, 2))data = np.vstack([cluster_1, cluster_2, cluster_3])# Aplicar K-Means con k=3 para obtener las etiquetas (usamos K-Means solo para obtener etiquetas iniciales)kmeans = KMeans(n_clusters=3, random_state=42).fit(data)labels = kmeans.labels_# Calcular los medoides manualmentemedoids = []for i in np.unique(labels): cluster_points = data[labels == i] distances = pairwise_distances(cluster_points, metric='euclidean') medoid_index = np.argmin(distances.sum(axis=0)) medoids.append(cluster_points[medoid_index])medoids = np.array(medoids)# Crear el gráficoplt.figure(figsize=(8, 6))sns.scatterplot(x=data[:, 0], y=data[:, 1], hue=labels, palette="viridis", style=labels, legend=None)plt.scatter(medoids[:, 0], medoids[:, 1], color='red', s=200, marker='*', label='Medoid')# Añadir las elipses alrededor de cada cluster (simplificadas)for i in np.unique(labels): cluster_points = data[labels == i] mean_x, mean_y = np.mean(cluster_points, axis=0) cov_matrix = np.cov(cluster_points.T) width =2* np.std(cluster_points[:, 0]) height =2* np.std(cluster_points[:, 1]) ellipse = Ellipse((mean_x, mean_y), width, height, edgecolor='black', facecolor='none', linewidth=2, linestyle='--') plt.gca().add_patch(ellipse)plt.title("Clusters con Medoides Destacados y Elipses Simplificadas")plt.xlabel("X")plt.ylabel("Y")plt.legend()# Guardar el gráfico en un archivo PNGplt.savefig("_figs/clusters_kmedoids.png", format="png", dpi=300, bbox_inches='tight')# Mostrar el gráficoplt.show()
Ejemplo LCS
Tiempo que demora esta sección: 0.1 minutos
Código
set.seed(324) # Setting seed for reproducibility # Slice specific rows, random selection with seed, and select columns "0" to "19" result_df_trim <- ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>%slice_sample(n =4) %>%# Randomly select 4 rows select("0":"19") # Select columns "0" to "19" (assuming they are in positions 1 to 20)
invisible("Me faltan 4: los de PAM con comb y seq para om y lcs (om2 y lcs2)")cqi_quarter<-rbind.data.frame(cbind.data.frame(algo="hac", type="om", time="quarter", k=2:15, corr=F, om_dist_quarter_c$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="hac", type="lcs", time="quarter", k=2:15, corr=F, lcs_dist_quarter_c$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="om", time="quarter", k=2:15, corr=F, pamRange_quarter_om$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="om", time="quarter", k=2:15, corr=T, pamRange_quarter_om2$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="lcs", time="quarter", k=2:15, corr=F, pamRange_quarter_lcs$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="lcs", time="quarter", k=2:15, corr=T, pamRange_quarter_lcs2$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))))|> dplyr::select(algo, type, time, k, corr, PBC, ASW, HC, HG, R2, R2sq)# round(summary(silhouette(as.integer(om_dist_quarter_c$clustering$cluster2), as.dist(dist_quarter_om)))$clus.avg.widths,2)[attr(rev(sort(table(om_dist_quarter_c$clustering$cluster2))),"names")]#functión para generar el gráficotabs_quarter_clus_sol<-rbind.data.frame(func_tab_range_clus(om_dist_quarter_c),func_tab_range_clus(lcs_dist_quarter_c),func_tab_range_clus(pamRange_quarter_om),func_tab_range_clus(pamRange_quarter_om2),func_tab_range_clus(pamRange_quarter_lcs),func_tab_range_clus(pamRange_quarter_lcs2))#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:##:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:# Inicializamos una lista para almacenar los resultadosresultados_list <-list()# Definimos un rango para los clusters a evaluarcluster_range <-2:15# Definimos los métodos y sus variablesmetodos <-list(hac_om =list(data = om_dist_quarter_c, dist = dist_quarter_om),hac_lcs =list(data = lcs_dist_quarter_c, dist = dist_quarter_lcs),pam_om0 =list(data = pamRange_quarter_om, dist = dist_quarter_om),pam_om1 =list(data = pamRange_quarter_om2, dist = dist_quarter_om),pam_lcs0 =list(data = pamRange_quarter_lcs, dist = dist_quarter_lcs),pam_lcs1 =list(data = pamRange_quarter_lcs2, dist = dist_quarter_lcs))# Número máximo de clusters para definir las columnasmax_clusters <-max(cluster_range)# Iteramos sobre cada métodofor (metodo innames(metodos)) {# Creamos un data frame temporal para cada método metodo_result <-data.frame()# Iteramos sobre cada cluster en el rangofor (cluster in cluster_range) {# Construimos el nombre del cluster dinámicamente cluster_name <-paste0("cluster", cluster)# Intentamos calcular los valores de silhouette silhouette_values <-tryCatch(round(summary(silhouette(as.integer(metodos[[metodo]]$data$clustering[[cluster_name]]), as.dist(metodos[[metodo]]$dist)))$clus.avg.widths[attr(rev(sort(table(metodos[[metodo]]$data$clustering[[cluster_name]]))),"names")], 2),error =function(e) rep(NA, cluster) )# Creamos un vector con las columnas llenando con NA si faltan valores silhouette_full <-c(silhouette_values, rep(NA, max_clusters -length(silhouette_values)))# Creamos un data frame temporal con los resultados para este cluster cluster_result <-data.frame(Metodo = metodo,Cluster = cluster,t(silhouette_full) # Transponemos los valores para que cada uno sea una columna )# Nombramos dinámicamente las columnas de silhouettecolnames(cluster_result)[3:(3+ max_clusters -1)] <-paste0("asw", 1:max_clusters)# Añadimos el resultado del cluster al data frame del método metodo_result <-rbind(metodo_result, cluster_result) }# Agregamos los resultados del método a la lista general resultados_list[[metodo]] <- metodo_result}# Combinamos todos los resultados en un único data frameavs_por_cluster_quarter <-do.call(rbind, resultados_list)# Ordenamos por Método y Clusteravs_por_cluster_quarter <- avs_por_cluster_quarter[order(avs_por_cluster_quarter$Metodo, avs_por_cluster_quarter$Cluster), ]bind_cols(cqi_quarter, tabs_quarter_clus_sol)%>% dplyr::mutate(corr= dplyr::case_when(corr==TRUE& algo!="hac"~"1",corr==FALSE& algo!="hac"~"0",T~""), key=paste0(algo,"_",type,corr,"_",k))|>left_join(dplyr::mutate(avs_por_cluster_quarter, key=paste0(Metodo,"_",Cluster)), by="key") |> dplyr::select(-Metodo, -Cluster) |>`rownames<-`(NULL) %>% dplyr::mutate(calc=round(PBC*(1/HC)*HG,2)) %>% dplyr::arrange(desc(ASW)) %>% dplyr::select(c("algo", "type", "time", "k", "corr", "PBC", "ASW", "HC", "HG", "R2", "R2sq", "calc", paste0("X",1:15), paste0("asw",1:15))) %>% {assign("asw_quarter_qci", dplyr::select(., -"time"), envir = .GlobalEnv) . } %>% { rio::export(., "_output/sol_conglomerados_tab_quarter.xlsx") knitr::kable(., "markdown", caption ="CQIs y frecuencias en conglomerados (trimestre)") }
CQIs y frecuencias en conglomerados (trimestre)
algo
type
time
k
corr
PBC
ASW
HC
HG
R2
R2sq
calc
X1
X2
X3
X4
X5
X6
X7
X8
X9
X10
X11
X12
X13
X14
X15
asw1
asw2
asw3
asw4
asw5
asw6
asw7
asw8
asw9
asw10
asw11
asw12
asw13
asw14
asw15
hac
lcs
quarter
2
0.61
0.70
0.11
0.94
0.10
0.23
5.21
5737
301
0.74
-0.06
pam
om
quarter
7
0
0.60
0.61
0.09
0.79
0.41
0.52
5.27
4372
680
313
205
189
142
137
0.70
0.64
0.21
0.52
0.06
-0.08
0.14
pam
om
quarter
13
0
0.61
0.61
0.06
0.86
0.50
0.59
8.74
3909
680
232
203
151
136
130
113
110
99
95
93
87
0.75
0.61
0.27
0.50
0.16
-0.13
0.08
0.33
-0.05
0.33
0.31
0.26
0.32
pam
om
quarter
14
0
0.61
0.61
0.05
0.87
0.51
0.60
10.61
3837
680
232
203
151
136
130
113
108
97
92
92
87
80
0.76
0.60
0.25
0.49
0.15
-0.14
0.07
0.32
-0.04
0.35
0.33
0.27
0.32
0.30
pam
om
quarter
15
0
0.61
0.61
0.05
0.86
0.52
0.60
10.49
3781
680
232
203
151
136
130
113
106
95
87
87
85
81
71
0.76
0.60
0.24
0.49
0.14
-0.14
0.06
0.31
-0.04
0.36
0.32
0.42
0.38
0.30
0.50
pam
om
quarter
4
0
0.50
0.60
0.18
0.67
0.31
0.39
1.86
4788
684
354
212
0.64
0.64
0.09
0.49
pam
om
quarter
5
0
0.57
0.60
0.13
0.74
0.36
0.47
3.24
4674
680
321
207
156
0.64
0.65
0.23
0.51
-0.05
pam
om
quarter
6
0
0.59
0.60
0.10
0.78
0.39
0.50
4.60
4476
680
319
207
202
154
0.68
0.64
0.22
0.50
0.02
-0.06
pam
om
quarter
12
0
0.61
0.60
0.07
0.85
0.49
0.58
7.41
3988
680
232
204
152
136
130
113
111
106
99
87
0.72
0.61
0.28
0.50
0.16
-0.13
0.09
0.35
-0.04
0.25
0.36
0.32
pam
om
quarter
4
1
0.50
0.60
0.18
0.67
0.31
0.39
1.86
4788
684
354
212
0.64
0.64
0.09
0.49
pam
om
quarter
5
1
0.57
0.60
0.13
0.74
0.36
0.47
3.24
4674
680
319
207
158
0.64
0.65
0.24
0.51
-0.05
pam
om
quarter
6
1
0.59
0.60
0.10
0.78
0.39
0.50
4.60
4476
680
318
207
203
154
0.68
0.64
0.23
0.50
0.01
-0.06
pam
om
quarter
7
1
0.61
0.60
0.08
0.82
0.41
0.52
6.25
4338
680
302
209
195
161
153
0.71
0.64
0.22
0.47
0.01
-0.04
-0.08
pam
om
quarter
12
1
0.61
0.60
0.07
0.85
0.49
0.58
7.41
3988
680
232
204
135
134
131
118
113
112
104
87
0.72
0.61
0.28
0.50
-0.13
0.16
0.08
0.46
0.22
-0.04
0.26
0.32
pam
om
quarter
15
1
0.61
0.60
0.05
0.87
0.53
0.61
10.61
3785
680
217
205
137
131
130
127
111
107
92
88
85
81
62
0.75
0.60
0.28
0.47
-0.06
0.13
0.04
0.37
0.17
0.05
0.61
0.20
0.33
0.54
-0.10
hac
om
quarter
4
0.49
0.59
0.20
0.65
0.29
0.37
1.59
4910
678
241
209
0.62
0.68
-0.10
0.53
hac
lcs
quarter
15
0.55
0.59
0.03
0.93
0.55
0.65
17.05
3117
1257
458
212
200
194
164
135
110
83
64
20
13
6
5
1.00
-0.16
1.00
0.47
-0.18
0.40
-0.03
0.17
-0.02
0.00
-0.05
0.03
0.05
0.47
-0.03
pam
om
quarter
8
0
0.59
0.59
0.09
0.80
0.43
0.53
5.24
4269
680
313
205
157
142
137
135
0.68
0.63
0.20
0.52
0.18
-0.08
0.13
0.24
pam
om
quarter
9
0
0.60
0.59
0.09
0.81
0.45
0.54
5.40
4164
680
313
204
156
142
137
132
110
0.70
0.63
0.18
0.52
0.17
-0.09
0.11
0.25
0.21
pam
om
quarter
10
0
0.61
0.59
0.08
0.84
0.46
0.56
6.40
4063
680
313
204
154
142
137
128
111
106
0.72
0.62
0.17
0.51
0.16
-0.11
0.10
0.26
-0.03
0.26
pam
om
quarter
11
0
0.61
0.59
0.07
0.84
0.48
0.57
7.32
3988
680
313
204
152
142
130
113
111
106
99
0.72
0.61
0.15
0.50
0.16
-0.12
0.11
0.35
-0.04
0.25
0.36
pam
om
quarter
8
1
0.61
0.59
0.08
0.82
0.43
0.53
6.25
4228
680
282
209
182
165
152
140
0.70
0.63
0.26
0.47
-0.08
0.04
-0.08
0.29
pam
om
quarter
10
1
0.61
0.59
0.08
0.84
0.47
0.56
6.40
4063
680
232
219
204
144
141
125
118
112
0.72
0.62
0.29
0.01
0.51
0.15
-0.11
0.40
0.15
-0.03
pam
om
quarter
11
1
0.61
0.59
0.07
0.84
0.48
0.57
7.32
3988
680
313
204
141
134
131
118
113
112
104
0.72
0.61
0.15
0.50
-0.12
0.16
0.09
0.46
0.22
-0.04
0.26
pam
om
quarter
14
1
0.61
0.59
0.05
0.87
0.51
0.60
10.61
3837
680
232
203
135
131
131
117
116
106
89
89
87
85
0.73
0.60
0.25
0.49
-0.14
0.35
0.06
0.02
0.17
-0.03
0.67
0.19
0.32
0.34
pam
om
quarter
9
1
0.60
0.58
0.09
0.81
0.45
0.54
5.40
4164
680
313
204
162
142
130
128
115
0.69
0.63
0.18
0.51
0.07
-0.09
0.13
0.39
0.18
pam
om
quarter
13
1
0.61
0.58
0.06
0.86
0.50
0.59
8.74
3909
680
232
203
135
131
130
114
114
109
105
89
87
0.71
0.61
0.26
0.50
-0.13
0.07
0.39
0.20
0.08
-0.05
0.08
0.67
0.32
hac
om
quarter
6
0.53
0.57
0.14
0.72
0.35
0.44
2.73
4502
678
408
226
209
15
0.64
0.66
-0.08
-0.07
0.51
0.22
hac
lcs
quarter
13
0.55
0.57
0.03
0.91
0.53
0.63
16.68
3117
1469
458
200
194
164
135
110
83
64
25
13
6
1.00
-0.15
1.00
-0.18
0.40
-0.03
0.41
-0.01
0.02
-0.05
-0.11
0.05
0.47
hac
lcs
quarter
14
0.55
0.57
0.03
0.91
0.54
0.64
16.68
3117
1469
458
200
194
164
135
110
83
64
20
13
6
5
1.00
-0.15
1.00
-0.18
0.40
-0.03
0.41
-0.01
0.02
-0.05
0.03
0.05
0.47
-0.03
pam
om
quarter
3
0
0.45
0.57
0.23
0.64
0.25
0.31
1.25
4982
688
368
0.59
0.64
0.07
pam
om
quarter
3
1
0.45
0.57
0.23
0.64
0.25
0.31
1.25
4982
688
368
0.59
0.64
0.07
hac
om
quarter
2
0.38
0.56
0.29
0.57
0.20
0.23
0.75
5151
887
0.58
0.44
hac
om
quarter
3
0.48
0.56
0.20
0.65
0.24
0.31
1.56
4910
887
241
0.62
0.41
-0.10
hac
om
quarter
5
0.53
0.56
0.14
0.72
0.34
0.41
2.73
4502
678
408
241
209
0.65
0.66
-0.07
-0.14
0.51
hac
om
quarter
15
0.58
0.56
0.08
0.82
0.49
0.60
5.94
3886
678
209
184
175
171
152
145
141
100
85
70
16
15
11
0.69
0.62
0.47
-0.05
0.53
0.13
0.23
-0.17
0.07
0.57
-0.09
-0.12
-0.15
0.10
0.06
pam
om
quarter
2
0
0.34
0.56
0.33
0.55
0.19
0.22
0.57
5348
690
0.55
0.65
pam
om
quarter
2
1
0.34
0.56
0.33
0.55
0.19
0.22
0.57
5348
690
0.55
0.65
pam
lcs
quarter
11
1
0.59
0.56
0.04
0.93
0.46
0.60
13.72
4098
638
269
192
179
168
153
134
98
66
43
0.67
0.57
0.32
-0.12
0.49
0.07
0.26
0.02
0.40
0.11
0.00
pam
lcs
quarter
12
1
0.59
0.56
0.03
0.95
0.47
0.61
18.68
4019
634
269
186
185
174
155
128
100
79
66
43
0.69
0.58
0.30
0.45
-0.21
0.04
0.28
0.03
0.46
0.36
0.11
0.00
pam
lcs
quarter
15
1
0.59
0.56
0.02
0.96
0.51
0.63
28.32
3817
638
269
192
177
168
134
109
108
102
80
69
66
66
43
0.70
0.54
0.28
-0.19
0.47
0.05
0.00
0.33
0.11
0.67
0.29
0.75
0.49
0.10
0.00
hac
om
quarter
9
0.58
0.55
0.10
0.78
0.40
0.52
4.52
4341
678
327
226
209
145
81
16
15
0.64
0.65
0.16
-0.11
0.50
-0.11
-0.16
-0.03
0.19
hac
om
quarter
10
0.58
0.55
0.10
0.80
0.42
0.54
4.64
4157
678
327
226
209
184
145
81
16
15
0.67
0.64
0.14
-0.13
0.49
0.03
-0.13
-0.17
-0.07
0.19
hac
om
quarter
11
0.58
0.55
0.10
0.80
0.43
0.55
4.64
4157
678
327
209
184
145
141
85
81
16
15
0.67
0.64
0.11
0.48
0.03
-0.13
0.12
-0.08
-0.17
-0.07
0.15
hac
om
quarter
12
0.58
0.55
0.10
0.80
0.45
0.57
4.64
4157
678
209
184
175
152
145
141
85
81
16
15
0.64
0.64
0.48
0.02
0.55
0.25
-0.15
0.10
-0.09
-0.20
-0.07
0.14
hac
lcs
quarter
11
0.54
0.55
0.04
0.91
0.51
0.61
12.29
3117
1469
658
194
164
135
110
96
64
25
6
1.00
-0.15
0.48
0.40
-0.03
0.41
0.00
-0.06
-0.04
-0.10
0.52
hac
lcs
quarter
12
0.54
0.55
0.04
0.91
0.51
0.62
12.29
3117
1469
658
194
164
135
110
83
64
25
13
6
1.00
-0.15
0.48
0.40
-0.03
0.41
-0.01
0.02
-0.05
-0.10
0.05
0.47
pam
lcs
quarter
10
1
0.58
0.55
0.05
0.92
0.44
0.56
10.67
4126
641
260
252
179
162
132
109
100
77
0.66
0.57
-0.15
0.39
0.49
0.07
-0.08
0.32
0.48
0.02
pam
lcs
quarter
13
1
0.59
0.55
0.03
0.94
0.49
0.61
18.49
3949
648
269
186
176
144
127
127
101
100
99
69
43
0.67
0.52
0.30
0.45
0.03
-0.16
0.03
0.09
0.67
0.34
0.50
0.08
0.00
hac
om
quarter
7
0.54
0.54
0.14
0.72
0.37
0.48
2.78
4502
678
327
226
209
81
15
0.59
0.66
0.19
-0.09
0.51
-0.14
0.19
hac
om
quarter
8
0.56
0.54
0.13
0.74
0.38
0.50
3.19
4486
678
327
226
209
81
16
15
0.60
0.66
0.18
-0.10
0.51
-0.14
-0.02
0.19
hac
om
quarter
13
0.58
0.54
0.08
0.82
0.47
0.58
5.94
3886
678
271
209
184
175
152
145
141
85
81
16
15
0.69
0.63
-0.05
0.47
-0.01
0.54
0.24
-0.17
0.08
-0.09
-0.21
-0.12
0.14
hac
om
quarter
14
0.58
0.54
0.08
0.82
0.47
0.59
5.94
3886
678
271
209
184
175
152
145
141
85
70
16
15
11
0.69
0.63
-0.05
0.47
-0.01
0.54
0.24
-0.17
0.08
-0.09
-0.12
-0.12
0.10
0.06
hac
lcs
quarter
9
0.53
0.54
0.05
0.87
0.49
0.58
9.22
3117
1633
658
194
174
135
96
25
6
1.00
-0.17
0.48
0.40
-0.09
0.42
-0.05
-0.01
0.59
hac
lcs
quarter
10
0.53
0.54
0.05
0.88
0.49
0.59
9.33
3117
1633
658
194
135
110
96
64
25
6
1.00
-0.17
0.48
0.40
0.42
0.00
-0.06
-0.04
-0.01
0.52
pam
lcs
quarter
6
1
0.54
0.54
0.12
0.84
0.36
0.47
3.78
4430
672
315
233
196
192
0.62
0.52
0.22
-0.06
0.43
0.11
pam
lcs
quarter
11
0
0.59
0.53
0.04
0.93
0.46
0.57
13.72
4046
675
253
206
181
179
132
110
99
81
76
0.65
0.47
-0.06
-0.05
0.62
0.48
-0.08
0.38
0.37
0.37
0.03
pam
lcs
quarter
7
1
0.56
0.53
0.09
0.87
0.39
0.51
5.41
4370
641
264
254
179
178
152
0.60
0.60
-0.09
0.40
0.04
0.53
-0.07
hac
lcs
quarter
8
0.52
0.52
0.07
0.84
0.44
0.55
6.24
3311
1633
658
174
135
96
25
6
0.89
-0.13
0.49
-0.09
0.42
-0.04
0.00
0.59
pam
lcs
quarter
6
0
0.54
0.52
0.12
0.84
0.36
0.47
3.78
4430
678
267
237
233
193
0.59
0.51
0.41
-0.06
0.00
0.45
pam
lcs
quarter
10
0
0.57
0.52
0.05
0.91
0.44
0.55
10.37
4077
677
274
206
181
178
152
110
101
82
0.64
0.47
-0.09
-0.04
0.62
0.50
-0.08
0.38
0.34
0.36
pam
lcs
quarter
15
0
0.58
0.52
0.02
0.95
0.51
0.63
27.55
3817
673
259
206
179
178
128
110
98
80
70
69
68
64
39
0.65
0.44
-0.08
-0.17
0.45
0.61
0.03
0.24
0.27
0.30
1.00
0.06
0.74
0.54
0.04
pam
lcs
quarter
4
1
0.54
0.52
0.15
0.82
0.28
0.41
2.95
4883
669
307
179
0.56
0.56
0.12
-0.05
hac
lcs
quarter
6
0.50
0.51
0.09
0.81
0.41
0.49
4.50
3311
1768
658
270
25
6
0.90
-0.13
0.49
-0.08
0.01
0.60
hac
lcs
quarter
7
0.51
0.51
0.07
0.84
0.43
0.52
6.12
3311
1633
658
270
135
25
6
0.89
-0.13
0.49
-0.15
0.43
0.00
0.60
pam
lcs
quarter
5
0
0.52
0.51
0.16
0.81
0.32
0.44
2.63
4611
678
314
218
217
0.57
0.52
0.25
0.07
-0.04
pam
lcs
quarter
9
0
0.57
0.51
0.06
0.90
0.42
0.53
8.55
4159
677
274
206
181
178
152
110
101
0.61
0.48
-0.08
-0.03
0.63
0.51
-0.07
0.40
0.37
pam
lcs
quarter
14
0
0.58
0.51
0.03
0.94
0.50
0.59
18.17
3841
675
253
206
181
179
132
110
99
81
76
70
69
66
0.63
0.44
-0.07
-0.17
0.60
0.45
-0.08
0.24
0.27
0.31
0.03
1.00
0.73
0.52
pam
lcs
quarter
5
1
0.52
0.51
0.16
0.81
0.32
0.43
2.63
4611
672
315
240
200
0.57
0.54
0.25
-0.07
0.12
pam
lcs
quarter
9
1
0.58
0.51
0.06
0.91
0.42
0.55
8.80
4223
621
260
222
177
172
152
132
79
0.59
0.62
-0.12
-0.01
0.66
0.53
0.16
-0.06
0.03
pam
lcs
quarter
14
1
0.59
0.51
0.02
0.95
0.50
0.62
28.02
3878
638
269
245
179
167
125
103
90
86
80
70
65
43
0.63
0.55
0.28
-0.23
0.46
0.05
0.04
0.31
0.43
0.26
0.40
1.00
0.11
-0.01
hac
lcs
quarter
4
0.49
0.50
0.09
0.80
0.39
0.43
4.36
3311
1768
658
301
0.90
-0.13
0.49
-0.13
hac
lcs
quarter
5
0.49
0.50
0.09
0.80
0.40
0.46
4.36
3311
1768
658
276
25
0.90
-0.13
0.49
-0.09
0.01
pam
lcs
quarter
3
0
0.46
0.50
0.21
0.74
0.23
0.32
1.62
5032
678
328
0.52
0.55
0.03
pam
lcs
quarter
4
0
0.50
0.50
0.18
0.79
0.27
0.37
2.19
4826
678
328
206
0.55
0.53
0.00
0.11
pam
lcs
quarter
8
0
0.56
0.50
0.08
0.88
0.40
0.52
6.16
4261
677
278
206
181
178
147
110
0.58
0.49
-0.08
0.00
0.64
0.51
-0.07
0.43
pam
lcs
quarter
13
0
0.58
0.50
0.03
0.94
0.48
0.59
18.17
3903
675
253
206
181
179
132
110
87
85
81
76
70
0.61
0.45
-0.06
-0.16
0.60
0.46
-0.08
0.25
0.31
0.55
0.32
0.03
1.00
pam
lcs
quarter
3
1
0.46
0.50
0.21
0.74
0.23
0.32
1.62
5032
678
328
0.52
0.55
0.03
pam
lcs
quarter
8
1
0.57
0.50
0.07
0.89
0.41
0.52
7.25
4259
636
264
251
187
178
149
114
0.58
0.60
-0.12
-0.05
0.45
0.65
-0.07
0.31
pam
lcs
quarter
7
0
0.56
0.48
0.09
0.87
0.38
0.50
5.41
4371
677
278
206
181
178
147
0.55
0.51
-0.07
0.04
0.65
0.53
-0.07
pam
lcs
quarter
12
0
0.58
0.48
0.04
0.93
0.47
0.58
13.48
3976
675
253
206
181
179
132
110
99
81
76
70
0.58
0.46
-0.06
-0.15
0.61
0.47
-0.08
0.25
0.36
0.36
0.03
1.00
hac
lcs
quarter
3
0.48
0.46
0.15
0.76
0.27
0.35
2.43
3969
1768
301
0.71
-0.02
-0.12
pam
lcs
quarter
2
0
0.22
0.45
0.43
0.51
0.13
0.14
0.26
5357
681
0.43
0.57
pam
lcs
quarter
2
1
0.22
0.45
0.43
0.51
0.13
0.14
0.26
5357
681
0.43
0.57
Tiempo que demora esta sección: 0 minutos
Se identificó una solución de 2, 4 y 7 conglomerados mediante el algoritmo de partición alrededor de medoides (PAM), utilizando distancias de emparejamiento óptimo (OM). Para validar la robustez de esta tipología, se implementó un procedimiento de bootstrap paramétrico con 1000 réplicas, comparando la calidad de la solución observada con la obtenida al aplicar el mismo procedimiento de clustering a datos generados bajo un modelo nulo combinado. Este modelo nulo evalúa la calidad de la agrupación en ausencia de estructura real, considerando aspectos combinados de duración y secuencia (comb), así como la secuencia por sí sola (seq), según la metodología propuesta por Studer (2021). Se busca medir cuánta calidad obtenida por la tipología sobrepasa la que habría sido obtenida para datos sin una estructura de conglomerados. El criterio Max T obtiene los máximos puntajes para las trayectorias que asumen estructuras aleatorias sobre las que comparar la solución de conglomerados obtenida.@studer_validating_2021
Se descartó la solución de 2 conglomerados por ser demasiado simple y, por lo tanto, carecer de valor explicativo en relación con los estados discretos de interés (a saber, las causas de hospitalización). En contraste con una estructura de trayectoria de duración y secuencia aleatoria, los valores obtenidos en los índices Gamma de Hubbert (HG) y Correlación Punto-Biserial (PBC) se encuentran dentro del intervalo esperado para la solución de 7 conglomerados, mientras que para la solución de 4 conglomerados, el índice C de Hubbert también se encuentra dentro de ese intervalo sumado a los otros. Esto significa que sólo el ASW de la solución de 4 conglomerados es mayor que lo esperado. Ahora, en comparación a una estructura de trayectorias aleatoria en términos secuencias, ambas soluciones tienen valores de ASW, HC y PBC superiores al esperable, mientras que en el caso del índice HG sólo es superior al esperable para la solución de 7 conglomerados.
A continuación se muestra con más detalle el resultado de pruebas de validación mediante bootstraps.
Código
# results: hide# fig.show: hideratio_plot=5asw_grob <-save_base_plot_as_grob(pam_om_quarter_null_comb_plot_asw, width =800*ratio_plot, height =600*ratio_plot, res=500)hc_grob <-save_base_plot_as_grob(pam_om_quarter_null_comb_plot_hc, width =800*ratio_plot, height =600*ratio_plot, res=500)hg_grob <-save_base_plot_as_grob(pam_om_quarter_null_comb_plot_hg, width =800*ratio_plot, height =600*ratio_plot, res=500)pbc_grob <-save_base_plot_as_grob(pam_om_quarter_null_comb_plot_pbc, width =800*ratio_plot, height =600*ratio_plot, res=500)final_plot_comb <-plot_grid( asw_grob, hc_grob, hg_grob, pbc_grob,ncol =2, # Número de columnasnrow =2, # Número de filasrel_widths =c(1, 1), # Ancho relativo de los gráficosrel_heights =c(1, 1),labels =c("A", "B", "C", "D"), # Etiquetas opcionaleslabel_size =15, # Tamaño de las etiquetasalign ="v", # Alineación vertical de los gráficosaxis ="tb"# Alineación de ejes superior e inferior)ggdraw() +draw_plot(final_plot_comb, x =0, y =0.1, width =1, height =0.9) +draw_text("Área gris: índices de agrupaciones aleatorias; línea negra: índices obtenidos", x =0.05, y =0.05, hjust =0, size =8, lineheight = .8)
Indicadores de calidad vs. bootstrap con secuencias y duraciones aleatorias
cqi_month<-rbind.data.frame(cbind.data.frame(algo="hac", type="om", time="month", k=2:15, corr=F, om_dist_month_c$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="hac", type="lcs", time="month", k=2:15, corr=F, lcs_dist_month_c$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="om", time="month", k=2:15, corr=F, pamRange_month_om$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="om", time="month", k=2:15, corr=T, pamRange_month_om2$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="lcs", time="month", k=2:15, corr=F, pamRange_month_lcs$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2))),cbind.data.frame(algo="pam", type="lcs", time="month", k=2:15, corr=T, pamRange_month_lcs2$stats) |> dplyr::mutate(across(PBC:HC,~round(.,2)))) |> dplyr::select(algo, type, time, k, corr, PBC, ASW, HC, HG, R2, R2sq)tabs_month_clus_sol<-rbind.data.frame(func_tab_range_clus(om_dist_month_c),func_tab_range_clus(lcs_dist_month_c),func_tab_range_clus(pamRange_month_om),func_tab_range_clus(pamRange_month_om2),func_tab_range_clus(pamRange_month_lcs),func_tab_range_clus(pamRange_month_lcs2))#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:#:# Inicializamos una lista para almacenar los resultadosresultados_list <-list()# Definimos un rango para los clusters a evaluarcluster_range <-2:15# Definimos los métodos y sus variablesmetodos <-list(hac_om =list(data = om_dist_month_c, dist = dist_month_om),hac_lcs =list(data = lcs_dist_month_c, dist = dist_month_lcs),pam_om0 =list(data = pamRange_month_om, dist = dist_month_om),pam_om1 =list(data = pamRange_month_om2, dist = dist_month_om),pam_lcs0 =list(data = pamRange_month_lcs, dist = dist_month_lcs),pam_lcs1 =list(data = pamRange_month_lcs2, dist = dist_month_lcs))# Número máximo de clusters para definir las columnasmax_clusters <-max(cluster_range)# Iteramos sobre cada métodofor (metodo innames(metodos)) {# Creamos un data frame temporal para cada método metodo_result <-data.frame()# Iteramos sobre cada cluster en el rangofor (cluster in cluster_range) {# Construimos el nombre del cluster dinámicamente cluster_name <-paste0("cluster", cluster)# Intentamos calcular los valores de silhouette silhouette_values <-tryCatch(round(summary(silhouette(as.integer(metodos[[metodo]]$data$clustering[[cluster_name]]), as.dist(metodos[[metodo]]$dist)))$clus.avg.widths[attr(rev(sort(table(metodos[[metodo]]$data$clustering[[cluster_name]]))),"names")], 2),error =function(e) rep(NA, cluster) )# Creamos un vector con las columnas llenando con NA si faltan valores silhouette_full <-c(silhouette_values, rep(NA, max_clusters -length(silhouette_values)))# Creamos un data frame temporal con los resultados para este cluster cluster_result <-data.frame(Metodo = metodo,Cluster = cluster,t(silhouette_full) # Transponemos los valores para que cada uno sea una columna )# Nombramos dinámicamente las columnas de silhouettecolnames(cluster_result)[3:(3+ max_clusters -1)] <-paste0("asw", 1:max_clusters)# Añadimos el resultado del cluster al data frame del método metodo_result <-rbind(metodo_result, cluster_result) }# Agregamos los resultados del método a la lista general resultados_list[[metodo]] <- metodo_result}# Combinamos todos los resultados en un único data frameavs_por_cluster_month <-do.call(rbind, resultados_list)# Ordenamos por Método y Clusteravs_por_cluster_month <- avs_por_cluster_month[order(avs_por_cluster_month$Metodo, avs_por_cluster_month$Cluster), ]bind_cols(cqi_month, tabs_month_clus_sol)%>% dplyr::mutate(corr= dplyr::case_when(corr==TRUE& algo!="hac"~"1",corr==FALSE& algo!="hac"~"0",T~""), key=paste0(algo,"_",type,corr,"_",k))|>left_join(dplyr::mutate(avs_por_cluster_month, key=paste0(Metodo,"_",Cluster)), by="key")|> dplyr::select(-Metodo, -Cluster) |>`rownames<-`(NULL) %>% dplyr::mutate(calc=round(PBC*(1/HC)*HG,2)) %>% dplyr::arrange(desc(ASW)) %>% dplyr::select(c("algo", "type", "time", "k", "corr", "PBC", "ASW", "HC", "HG", "R2", "R2sq", "calc", paste0("X",1:15), paste0("asw",1:15)))%>% {assign("asw_month_qci", dplyr::select(., -"time"), envir = .GlobalEnv) . } %>% { rio::export(., "_output/sol_conglomerados_tab_month.xlsx") knitr::kable(., "markdown", caption ="CQIs y frecuencias en conglomerados (mensual)") }
CQIs y frecuencias en conglomerados (mensual)
algo
type
time
k
corr
PBC
ASW
HC
HG
R2
R2sq
calc
X1
X2
X3
X4
X5
X6
X7
X8
X9
X10
X11
X12
X13
X14
X15
asw1
asw2
asw3
asw4
asw5
asw6
asw7
asw8
asw9
asw10
asw11
asw12
asw13
asw14
asw15
hac
om
month
2
0.51
0.85
0.17
0.98
0.04
0.13
2.94
5998
40
0.85
-0.06
hac
lcs
month
2
0.53
0.60
0.14
0.85
0.15
0.21
3.22
5174
864
0.74
-0.22
hac
lcs
month
3
0.56
0.56
0.11
0.86
0.18
0.32
4.38
5174
758
106
0.67
-0.09
-0.20
hac
lcs
month
4
0.56
0.56
0.11
0.86
0.19
0.38
4.38
5174
758
100
6
0.67
-0.09
-0.18
0.19
hac
om
month
9
0.43
0.48
0.18
0.66
0.29
0.44
1.58
4595
683
301
203
179
37
35
3
2
0.54
0.56
0.08
-0.20
0.38
-0.12
-0.22
0.52
0.75
hac
om
month
10
0.43
0.48
0.18
0.66
0.29
0.46
1.58
4595
683
301
203
179
35
29
8
3
2
0.54
0.56
0.08
-0.20
0.38
-0.22
-0.10
0.34
0.46
0.75
hac
om
month
11
0.43
0.48
0.18
0.66
0.30
0.48
1.58
4595
683
301
203
179
32
29
8
3
3
2
0.54
0.56
0.08
-0.20
0.38
-0.17
-0.10
0.34
0.44
0.23
0.75
pam
om
month
15
0
0.45
0.48
0.10
0.76
0.40
0.51
3.42
3889
681
542
180
142
97
87
85
80
69
65
62
28
26
5
0.59
0.51
0.35
0.22
0.23
0.11
0.01
-0.11
-0.18
-0.14
-0.01
-0.20
0.00
-0.13
0.21
pam
om
month
13
1
0.44
0.48
0.11
0.75
0.39
0.50
3.00
3909
682
555
203
181
107
87
80
76
71
65
16
6
0.60
0.51
0.31
0.11
0.22
-0.17
-0.06
-0.17
0.08
-0.19
0.02
0.04
0.15
pam
om
month
15
1
0.45
0.48
0.10
0.76
0.41
0.53
3.42
3889
681
542
199
180
94
81
74
71
64
61
59
25
13
5
0.59
0.51
0.35
0.17
0.21
-0.12
-0.18
0.04
0.04
-0.17
0.02
-0.11
-0.12
0.07
0.17
pam
om
month
9
0
0.42
0.47
0.14
0.70
0.36
0.45
2.10
4076
687
571
215
184
123
87
80
15
0.55
0.51
0.30
0.12
0.24
-0.19
-0.03
-0.18
0.01
pam
om
month
10
0
0.43
0.47
0.13
0.72
0.37
0.46
2.38
4016
684
559
183
159
134
108
103
77
15
0.57
0.51
0.31
0.23
0.16
-0.05
-0.18
-0.06
-0.19
0.02
pam
om
month
12
0
0.43
0.47
0.12
0.73
0.38
0.49
2.62
3995
683
538
183
142
124
103
96
74
64
30
6
0.56
0.51
0.36
0.22
0.20
0.10
-0.09
-0.15
-0.19
-0.20
-0.03
0.19
pam
om
month
13
0
0.44
0.47
0.11
0.75
0.39
0.49
3.00
3918
683
538
183
142
124
103
95
79
73
64
30
6
0.58
0.51
0.35
0.21
0.20
0.10
-0.10
-0.15
-0.17
-0.19
-0.20
-0.03
0.19
pam
om
month
14
0
0.44
0.47
0.10
0.76
0.39
0.50
3.34
3895
682
547
180
145
99
89
85
80
70
65
65
30
6
0.60
0.51
0.32
0.22
0.15
0.08
-0.04
-0.11
-0.17
-0.14
-0.22
0.02
-0.03
0.20
pam
om
month
9
1
0.42
0.47
0.14
0.71
0.35
0.43
2.13
4057
686
572
184
177
142
121
78
21
0.56
0.51
0.28
0.23
0.12
-0.06
-0.18
-0.20
-0.07
pam
om
month
10
1
0.43
0.47
0.13
0.72
0.37
0.46
2.38
4016
684
559
183
159
134
108
103
77
15
0.57
0.51
0.31
0.23
0.16
-0.05
-0.18
-0.06
-0.19
0.02
pam
om
month
11
1
0.43
0.47
0.13
0.72
0.37
0.48
2.38
4016
684
559
183
147
127
105
104
77
30
6
0.57
0.51
0.31
0.23
0.14
0.07
-0.13
-0.09
-0.21
-0.01
0.20
pam
om
month
12
1
0.43
0.47
0.12
0.73
0.38
0.50
2.62
3986
682
555
203
181
105
92
77
71
64
16
6
0.57
0.51
0.33
0.11
0.24
-0.17
-0.08
0.08
-0.19
0.04
0.04
0.15
pam
om
month
14
1
0.44
0.47
0.10
0.76
0.40
0.51
3.34
3895
682
547
203
180
91
90
79
67
63
61
60
14
6
0.60
0.51
0.32
0.09
0.22
-0.05
-0.12
-0.17
-0.19
0.07
0.06
-0.11
0.06
0.15
hac
om
month
8
0.43
0.46
0.18
0.66
0.26
0.42
1.58
4595
862
301
203
37
35
3
2
0.55
0.34
0.09
-0.19
-0.12
-0.22
0.52
0.77
pam
om
month
6
0
0.39
0.46
0.18
0.64
0.31
0.36
1.39
4190
688
674
208
190
88
0.53
0.51
0.25
-0.03
0.23
-0.24
pam
om
month
7
0
0.41
0.46
0.16
0.68
0.33
0.38
1.74
4105
685
658
198
188
120
84
0.55
0.50
0.28
-0.03
0.23
-0.11
-0.25
pam
om
month
8
0
0.42
0.46
0.14
0.70
0.34
0.38
2.10
4076
687
571
213
184
128
91
88
0.55
0.51
0.30
0.11
0.24
-0.21
-0.14
-0.26
pam
om
month
11
0
0.43
0.46
0.12
0.73
0.38
0.47
2.62
3998
683
545
183
157
133
102
95
77
59
6
0.56
0.51
0.35
0.23
0.18
-0.05
-0.06
-0.12
-0.19
-0.22
0.20
pam
om
month
7
1
0.40
0.46
0.17
0.65
0.33
0.37
1.53
4164
686
583
190
175
147
93
0.54
0.52
0.25
0.23
0.18
-0.08
-0.26
pam
om
month
8
1
0.42
0.46
0.14
0.70
0.34
0.38
2.10
4076
687
571
217
184
124
91
88
0.55
0.51
0.30
0.08
0.24
-0.20
-0.13
-0.26
hac
om
month
4
0.37
0.45
0.25
0.58
0.21
0.28
0.86
4798
864
336
40
0.51
0.35
-0.07
-0.10
hac
om
month
5
0.37
0.45
0.25
0.58
0.21
0.32
0.86
4798
864
336
37
3
0.51
0.35
-0.07
-0.08
0.52
hac
om
month
6
0.38
0.45
0.25
0.58
0.22
0.35
0.88
4798
862
336
37
3
2
0.50
0.37
-0.08
-0.09
0.52
0.77
hac
om
month
7
0.38
0.45
0.24
0.58
0.23
0.39
0.92
4798
862
301
37
35
3
2
0.49
0.37
0.11
-0.10
-0.21
0.52
0.77
pam
om
month
5
0
0.36
0.45
0.21
0.61
0.28
0.30
1.05
4225
698
693
230
192
0.53
0.19
0.49
-0.17
0.21
pam
om
month
5
1
0.36
0.45
0.21
0.61
0.28
0.29
1.05
4225
693
613
315
192
0.53
0.49
0.19
-0.14
0.21
pam
om
month
6
1
0.37
0.45
0.20
0.62
0.30
0.33
1.15
4215
693
610
230
192
98
0.54
0.48
0.20
0.01
0.21
-0.16
pam
lcs
month
2
0
0.26
0.45
0.38
0.50
0.10
0.13
0.34
5145
893
0.53
-0.05
pam
lcs
month
2
1
0.26
0.45
0.38
0.50
0.10
0.13
0.34
5145
893
0.53
-0.05
pam
lcs
month
14
0
0.43
0.44
0.09
0.82
0.41
0.54
3.92
3851
677
534
170
159
114
101
77
76
75
71
58
57
18
0.56
0.39
0.31
0.14
0.35
-0.04
0.04
-0.18
-0.14
0.18
-0.12
-0.14
0.07
-0.02
pam
lcs
month
11
1
0.41
0.44
0.12
0.78
0.38
0.51
2.66
3938
681
627
168
166
115
101
85
82
57
18
0.55
0.38
0.30
0.12
0.17
-0.02
-0.15
0.14
-0.15
-0.04
0.00
pam
lcs
month
15
1
0.43
0.44
0.09
0.83
0.41
0.55
3.97
3843
670
527
167
166
121
90
84
73
68
67
58
55
35
14
0.56
0.41
0.36
0.16
0.15
-0.08
0.04
-0.10
-0.14
0.23
-0.03
-0.15
-0.09
-0.14
0.02
hac
om
month
3
0.28
0.43
0.36
0.45
0.16
0.23
0.35
5134
864
40
0.44
0.39
-0.07
hac
om
month
14
0.41
0.43
0.14
0.68
0.36
0.52
1.99
3967
683
613
265
203
179
36
32
29
15
8
3
3
2
0.54
0.53
0.07
0.10
-0.24
0.35
0.16
-0.20
-0.12
-0.17
0.27
0.44
0.23
0.75
hac
om
month
15
0.41
0.43
0.14
0.68
0.36
0.53
1.99
3967
683
613
265
179
149
54
36
32
29
15
8
3
3
2
0.54
0.53
0.07
0.10
0.35
-0.15
-0.16
0.16
-0.20
-0.14
-0.17
0.26
0.44
0.23
0.75
pam
om
month
3
0
0.30
0.43
0.28
0.53
0.22
0.21
0.57
4440
900
698
0.52
-0.06
0.49
pam
om
month
4
0
0.33
0.43
0.25
0.56
0.25
0.26
0.74
4404
706
698
230
0.48
0.21
0.49
-0.16
pam
om
month
3
1
0.30
0.43
0.28
0.53
0.22
0.21
0.57
4440
900
698
0.52
-0.06
0.49
pam
lcs
month
11
0
0.41
0.43
0.12
0.78
0.38
0.51
2.66
3938
684
552
175
168
147
103
102
75
72
22
0.55
0.37
0.24
0.41
0.14
-0.12
-0.16
0.02
-0.13
0.02
-0.02
pam
lcs
month
12
0
0.42
0.43
0.10
0.79
0.39
0.53
3.32
3924
682
538
180
171
115
90
89
86
77
65
21
0.53
0.37
0.34
0.31
0.14
-0.04
0.12
-0.12
0.01
-0.13
-0.18
-0.03
pam
lcs
month
13
0
0.42
0.43
0.10
0.80
0.40
0.53
3.36
3890
684
539
171
159
116
101
79
77
75
71
58
18
0.55
0.36
0.30
0.12
0.35
-0.05
0.04
-0.14
-0.18
0.18
-0.12
0.06
-0.02
pam
lcs
month
15
0
0.43
0.43
0.08
0.83
0.42
0.55
4.46
3830
676
519
170
158
113
101
75
74
74
71
58
53
48
18
0.55
0.39
0.39
0.13
0.33
-0.05
0.04
0.17
-0.13
-0.20
-0.12
-0.15
0.11
-0.10
-0.01
pam
lcs
month
12
1
0.42
0.43
0.10
0.80
0.39
0.52
3.36
3903
680
541
190
163
124
110
83
83
79
61
21
0.55
0.38
0.33
0.22
0.18
-0.08
-0.06
-0.16
-0.18
0.18
0.04
-0.02
pam
lcs
month
13
1
0.43
0.43
0.10
0.80
0.40
0.55
3.44
3904
678
539
189
163
120
111
81
80
79
71
18
5
0.54
0.39
0.34
0.21
0.17
-0.07
-0.06
-0.15
0.18
-0.09
-0.19
0.09
0.23
pam
lcs
month
14
1
0.43
0.43
0.09
0.82
0.41
0.55
3.92
3869
670
525
165
163
131
94
88
80
76
63
60
38
16
0.54
0.42
0.38
0.22
0.17
-0.08
0.10
0.05
-0.09
-0.14
-0.18
-0.15
-0.01
0.00
hac
om
month
12
0.40
0.42
0.15
0.67
0.34
0.49
1.79
3982
683
613
301
203
179
32
29
8
3
3
2
0.53
0.53
0.09
0.00
-0.23
0.36
-0.19
-0.11
0.34
0.44
0.23
0.75
hac
om
month
13
0.40
0.42
0.15
0.67
0.35
0.50
1.79
3982
683
613
265
203
179
36
32
29
8
3
3
2
0.53
0.53
0.08
0.10
-0.24
0.36
0.16
-0.20
-0.12
0.27
0.44
0.23
0.75
hac
lcs
month
5
0.46
0.42
0.13
0.78
0.27
0.42
2.76
4507
758
667
100
6
0.51
-0.13
0.49
-0.18
0.19
hac
lcs
month
6
0.46
0.42
0.13
0.78
0.28
0.45
2.76
4507
758
667
88
12
6
0.51
-0.13
0.49
-0.16
0.15
0.13
hac
lcs
month
7
0.46
0.42
0.13
0.78
0.29
0.47
2.76
4507
758
667
88
10
6
2
0.51
-0.13
0.49
-0.17
0.39
0.09
0.76
pam
om
month
4
1
0.33
0.42
0.25
0.56
0.25
0.26
0.74
4404
698
621
315
0.48
0.49
0.20
-0.14
pam
lcs
month
7
0
0.39
0.42
0.16
0.72
0.33
0.42
1.75
4067
687
649
235
169
126
105
0.50
0.39
0.32
-0.03
0.17
-0.02
-0.19
pam
lcs
month
10
0
0.39
0.42
0.14
0.74
0.36
0.49
2.06
4007
685
548
170
163
146
111
101
85
22
0.51
0.38
0.23
0.15
-0.18
0.54
-0.02
-0.14
0.17
-0.01
pam
lcs
month
8
1
0.39
0.42
0.15
0.72
0.34
0.48
1.87
4054
682
647
235
170
119
109
22
0.51
0.40
0.30
-0.02
0.16
0.01
-0.14
0.00
pam
lcs
month
9
1
0.39
0.42
0.14
0.74
0.35
0.48
2.06
3994
679
641
226
165
111
101
99
22
0.52
0.39
0.31
-0.04
0.17
-0.01
-0.14
0.03
-0.01
pam
lcs
month
10
1
0.39
0.42
0.14
0.74
0.36
0.49
2.06
4009
682
554
166
154
130
122
103
98
20
0.51
0.40
0.22
0.17
0.51
-0.14
-0.05
-0.15
0.10
0.01
hac
lcs
month
15
0.37
0.41
0.06
0.81
0.45
0.61
5.00
2635
1709
667
278
269
163
134
77
36
27
25
10
3
3
2
1.00
-0.21
0.30
-0.20
0.22
0.15
-0.16
0.03
-0.09
-0.05
-0.12
0.27
0.45
0.27
0.75
pam
om
month
2
0
0.16
0.41
0.48
0.39
0.12
0.11
0.13
5337
701
0.40
0.54
pam
om
month
2
1
0.16
0.41
0.48
0.39
0.12
0.11
0.13
5337
701
0.40
0.54
pam
lcs
month
6
0
0.38
0.41
0.17
0.69
0.31
0.37
1.54
4115
687
685
238
173
140
0.50
0.40
0.23
-0.05
0.16
-0.19
pam
lcs
month
8
0
0.39
0.41
0.15
0.71
0.34
0.43
1.85
4078
687
563
198
172
125
119
96
0.49
0.39
0.28
0.33
0.17
-0.14
-0.01
-0.20
pam
lcs
month
9
0
0.39
0.41
0.15
0.72
0.35
0.48
1.87
4066
686
553
171
169
146
120
105
22
0.50
0.39
0.22
0.16
-0.16
0.55
0.01
-0.13
-0.01
hac
lcs
month
12
0.37
0.40
0.08
0.78
0.41
0.58
3.61
2798
1709
667
412
269
77
63
25
10
3
3
2
0.89
-0.16
0.33
-0.22
0.23
0.04
-0.14
-0.11
0.30
0.45
0.30
0.76
hac
lcs
month
13
0.37
0.40
0.08
0.78
0.42
0.59
3.61
2798
1709
667
412
269
77
36
27
25
10
3
3
2
0.89
-0.16
0.33
-0.22
0.23
0.04
-0.06
-0.04
-0.11
0.27
0.45
0.27
0.76
hac
lcs
month
14
0.37
0.40
0.07
0.79
0.43
0.60
4.18
2798
1709
667
278
269
134
77
36
27
25
10
3
3
2
0.89
-0.16
0.33
-0.19
0.23
-0.15
0.03
-0.09
-0.05
-0.12
0.27
0.45
0.27
0.76
hac
lcs
month
10
0.36
0.39
0.08
0.77
0.40
0.54
3.46
2798
1709
667
412
346
63
25
10
6
2
0.89
-0.16
0.33
-0.19
0.04
-0.12
-0.11
0.36
0.07
0.76
hac
lcs
month
11
0.36
0.39
0.08
0.77
0.40
0.56
3.46
2798
1709
667
412
346
63
25
10
3
3
2
0.89
-0.16
0.33
-0.19
0.04
-0.12
-0.11
0.34
0.45
0.30
0.76
pam
lcs
month
5
0
0.37
0.39
0.19
0.65
0.28
0.35
1.27
4288
687
685
238
140
0.45
0.40
0.25
-0.05
-0.19
hac
lcs
month
8
0.35
0.38
0.10
0.74
0.37
0.50
2.59
2798
1709
758
667
88
10
6
2
0.89
-0.16
-0.16
0.33
-0.17
0.39
0.09
0.76
hac
lcs
month
9
0.35
0.38
0.10
0.74
0.38
0.52
2.59
2798
1709
758
667
63
25
10
6
2
0.89
-0.16
-0.16
0.33
-0.10
-0.11
0.37
0.07
0.76
pam
lcs
month
4
0
0.33
0.38
0.23
0.60
0.24
0.28
0.86
4368
719
693
258
0.44
0.20
0.38
-0.19
pam
lcs
month
3
0
0.26
0.36
0.30
0.52
0.19
0.18
0.45
4451
893
694
0.44
-0.09
0.40
pam
lcs
month
3
1
0.26
0.36
0.30
0.52
0.19
0.18
0.45
4451
893
694
0.44
-0.09
0.40
pam
lcs
month
4
1
0.34
0.33
0.22
0.61
0.26
0.31
0.94
4380
687
624
347
0.34
0.42
0.47
-0.24
pam
lcs
month
5
1
0.37
0.32
0.19
0.66
0.29
0.35
1.29
4288
681
570
348
151
0.32
0.42
0.57
-0.16
-0.21
pam
lcs
month
7
1
0.39
0.32
0.15
0.72
0.34
0.42
1.87
4067
681
502
386
166
129
107
0.33
0.41
0.70
-0.17
0.17
-0.04
-0.19
pam
lcs
month
6
1
0.38
0.29
0.17
0.69
0.31
0.40
1.54
4233
681
502
386
129
107
0.28
0.41
0.71
-0.17
-0.03
-0.19
Código
# hac lcs month 3, pocos conglomerados, pero poco informativo# pam om month 6 FALSO
Tiempo que demora esta sección: 0 minutos
Código
frobenius_norm(as.matrix(func_tab_range_clus(pamRange_quarter_lcs)), as.matrix(func_tab_range_clus(pamRange_quarter_lcs2)))invisible("Frobenius norm: ||A-B||_F=\sqrt{\sum{i,j}(A_{ij}-B_ij)^2}")invisible("Por ahora lo dejamos pasar, es para comparar matrices y sus diferencias")
Tiempo que demora esta sección: 0 minutos
La mayoría de las soluciones presentaron conglomerados con tamaños muy pequeños (n < 30), lo que limita su generalización y sugiere una estabilidad presumiblemente baja, o bien con valores ASW negativos, lo que añade evidencia a problemas de estabilidad de la solución.
Código
asw_month_qci$label<-with(asw_month_qci, paste0(algo, "_",type, "_", k, "_", ifelse(corr==1,1,ifelse(corr==0,0,0))))plot(sort(asw_month_qci$ASW, decreasing =TRUE), type ="b", pch =19, xlab ="", xaxt ="n",ylab ="Valor de ASW", main ="", col ="blue")axis(1, at =1:length(asw_month_qci$label), labels = asw_month_qci$label, las =2, cex.axis =0.60)
Gráfico de codo, ASW con etiquetas soluciones de conglomerados
Tiempo que demora esta sección: 0 minutos
Se decidió por una solución de dos conglomerados que aunque tuvo valores globales en ASW subóptimos (0,41) y podría no ser del todo compleja, fue la única solución que obtuvo valores positvios en las medidas ASW por cada conglomerado, y sus conglomerados están integrados por más de 30 observaciones. A continuación, vemos los índices de calidad de la solución mediante remuestreos bootstrap con al menos 1,000 replicaciones.
Sólo los índices ASW y HC (C de Hubbert) mostraron mejores índices de calidad que los esperados para una estructura aleatoria en términos de secuencias.
A continuación se muestra el resultado de pruebas de validación mediante bootstraps.
Código
# results: hide# fig.show: hideratio_plot=5asw_grob_m <-save_base_plot_as_grob(pam_om_month_null_comb_plot_asw, width =800*ratio_plot, height =600*ratio_plot, res=500)hc_grob_m <-save_base_plot_as_grob(pam_om_month_null_comb_plot_hc, width =800*ratio_plot, height =600*ratio_plot, res=500)hg_grob_m <-save_base_plot_as_grob(pam_om_month_null_comb_plot_hg, width =800*ratio_plot, height =600*ratio_plot, res=500)pbc_grob_m <-save_base_plot_as_grob(pam_om_month_null_comb_plot_pbc, width =800*ratio_plot, height =600*ratio_plot, res=500)final_plot_comb_m <-plot_grid( asw_grob_m, hc_grob_m, hg_grob_m, pbc_grob_m,ncol =2, # Número de columnasnrow =2, # Número de filasrel_widths =c(1, 1), # Ancho relativo de los gráficosrel_heights =c(1, 1),labels =c("A", "B", "C", "D"), # Etiquetas opcionaleslabel_size =15, # Tamaño de las etiquetasalign ="v", # Alineación vertical de los gráficosaxis ="tb"# Alineación de ejes superior e inferior)ggdraw() +draw_plot(final_plot_comb_m, x =0, y =0.1, width =1, height =0.9) +draw_text("Área gris: índices de agrupaciones aleatorias; línea negra: índices obtenidos", x =0.05, y =0.05, hjust =0, size =8, lineheight = .8)
Indicadores de calidad vs. bootstrap con secuencias y duraciones aleatorias
invisible("Información sobre la solución")invisible("H:/Mi unidad/PERSONAL ANDRES/UCH_salud_publica/asignaturas/un_inv_II/_hist_sintaxis/un_inv_ii5_explorar_soluciones.R")# invisible("Hacemos clasificación de pertenencia cluster y ponemos etiquetas")ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4 <-factor(pamRange_quarter_om$clustering$cluster4,levels=rev(attr( sort(table(pamRange_quarter_om$clustering$cluster4)), "name")),labels=c("6035, Un trimestre, TSM(4)", "6025, Un trimestre, TUS(3)", "5939, Un semestre TSM(1)", "5989, Comorbilidad un trimestre(2)"))invisible("Me da buena: 0,61 en promedio. El problema está con 5710 y 6036 que son pésimos")sil_pam_om_clus4_q_nostd<-silhouette(as.integer(pamRange_quarter_om$clustering$cluster4), as.dist(dist_quarter_om))# Crear etiquetas personalizadascluster_labels <-paste0("Cluster ", seq_along(attr(summary(sil_pam_om_clus4_q_nostd)$clus.avg.widths, "dimnames")[[1]]), ":\nAWS ", sprintf("%1.2f",summary(sil_pam_om_clus4_q_nostd)$clus.avg.widths))# Graficar con etiquetas personalizadasfviz_silhouette( sil_pam_om_clus4_q_nostd, lab.clusters = cluster_labels, # Etiquetas personalizadas para los clústeresprint.summary=F) +scale_fill_grey(start =0.2, end =0.8, labels = cluster_labels) +# Escala de grisesscale_color_grey(start =0.2, end =0.8, labels = cluster_labels)+# Escala de grises para los bordeaggtitle(NULL)+labs(y="Ancho medio de la silueta", x="Conglomerados")# Elimina el título
length(subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, clus_pam_om4=="5939, Un semestre TSM(1)")$run)
[1] 354
Código
diag_pam_om4_5939<-df_filled %>% dplyr::filter(run %in%subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, clus_pam_om4=="5939, Un semestre TSM(1)")$run) %>% dplyr::select(run, diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11, fecha_egreso_rec_fmt, estab_homo) %>% dplyr::group_by(run) %>% dplyr::filter(row_number() !=1) %>%# Elimina la primera observación de cada run dplyr::mutate(all_diags =paste(na.omit(c(diag1, diag2, diag3, diag4, diag5, diag6, diag7, diag8, diag9, diag10, diag11)), collapse =", ") ) %>% dplyr::summarise(all_diags =first(all_diags),fecha_egreso_rec_fmt =first(fecha_egreso_rec_fmt),estab_homo =first(estab_homo) ) %>% dplyr::ungroup() %>% dplyr::pull(all_diags) %>%# Extraer la columna all_diags como vectorstrsplit(split =", ") %>%# Separar cada diagnóstico por comasunlist() # Convertir la lista en un vectorinvisible("head(arrange(data.frame(table(diag_pam_om4_5939)),-Freq),20)")invisible("Para chatgpt= estos son códigos de CIE-10, descríbeme brevemente cada uno en markdown en formato 'Cód. CIE-10 (n=Freq) - [descripción] '")cat("número de diagnósticos distintos")
número de diagnósticos distintos
Código
length(diag_pam_om4_5939)
[1] 1838
Código
# ggplot( head(arrange(data.frame(table(diag_pam_om4_5939)),-Freq),20), aes(x = reorder(diag_pam_om4_5939, -Freq), y = Freq)) +# geom_bar(stat = "identity", fill = "steelblue") +# labs(title = "Top 20 Most Frequent Diagnoses",# x = "Diagnosis",# y = "Frequency") +# theme(axis.text.x = element_text(angle = 45, hjust = 1))## F329 (n=131) - Episodio depresivo no especificado # F603 (n=91) - Trastorno de la personalidad emocionalmente inestable, tipo límite # F609 (n=91) - Otros trastornos específicos de la personalidad # F322 (n=86) - Episodio depresivo grave sin síntomas psicóticos # F209 (n=62) - Esquizofrenia no especificada # F192 (n=51) - Trastornos mentales y del comportamiento debidos al uso de múltiples drogas y al uso de otras sustancias psicoactivas # F200 (n=45) - Esquizofrenia paranoide # F319 (n=43) - Trastorno depresivo recurrente, episodio actual no especificado # F432 (n=32) - Trastornos de adaptación # Z915 (n=31) - Antecedentes personales de traumatismo no clasificado en otra parte # C490 (n=29) - Neoplasia maligna del tejido conectivo y de los tejidos blandos, de localización no especificada # G909 (n=29) - Trastorno del sistema nervioso no especificado # F323 (n=27) - Episodio depresivo grave con síntomas psicóticos # F431 (n=27) - Reacción al estrés agudo y trastornos de adaptación # F29X (n=21) - Psicosis no orgánica no especificada # Z511 (n=20) - Atención sanitaria para radioterapia # F419 (n=18) - Trastorno de ansiedad no especificado # F608 (n=16) - Otros trastornos específicos de la personalidad # F239 (n=14) - Trastorno psicótico agudo y transitorio no especificado # F449 (n=14) - Trastorno de ansiedad fóbica no especificado #
Tiempo que demora esta sección: 0 minutos
Entre quienes se encuentren en un semestre en el sistema por TSM y presentan un segundo episodio, de las 20 condiciones más frecuentes, al menos el 50% se caracteriza por episodios depresivos no especificado (F329), trastornos de la personalidad tipo límite (F603), otros trastornos específicos de la personalidad (F609) y episodios depresivos graves sin síntomas psicóticos (F322) y esquizofrenia no especificada (F209).
# X-squared = 136.81, df = 24, p-value < 2.2e-16invisible("No exclusivo de la clasificación por conglomerados, sino de las trayectorias en general y las matrices de distancia")#https://cran.r-project.org/web/packages/TraMineRextras/TraMineRextras.pdfseqcompare_sex_quarter_om<-seqCompare(States_Wide.seq_quarter_t_prim_adm, group=ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$glosa_sexo, method="OM", sm=costmatrix_quarter, seed=2125, s=1e1, stat="all")# LRT p-value BIC diff. Bayes Factor (Avg) Bayes Factor (From Avg BIC)# [1,] 3.723906 0.07729539 0.7281733 2.140152 1.439199seqcompare_sex_quarter_lcs<-seqCompare(States_Wide.seq_quarter_t_prim_adm, group=ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$glosa_sexo, method="LCS", sm=costmatrix_quarter, seed=2125, s=1e1, stat="all")# LRT p-value BIC diff. Bayes Factor (Avg) Bayes Factor (From Avg BIC)# [1,] 3.537305 0.08879071 0.5415731 2.120956 1.310995seqcompare_ppoo_quarter_om<-seqCompare(States_Wide.seq_quarter_t_prim_adm, group=ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$inclusivo_real_historico, method="OM", sm=costmatrix_quarter, seed=2125, s=1e1, stat="all")# Currently seqLRT supports only 2 groups!# LRT p-value BIC diff. Bayes Factor (Avg) Bayes Factor (From Avg BIC)# [1,] 3.321524 0.08714022 0.325792 1.479436 1.176914seqcompare_ppoo_quarter_lcs<-seqCompare(States_Wide.seq_quarter_t_prim_adm, group=ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$inclusivo_real_historico, method="LCS", sm=costmatrix_quarter, seed=2125, s=1e1, stat="all")# LRT p-value BIC diff. Bayes Factor (Avg) Bayes Factor (From Avg BIC)# [1,] 3.214252 0.09793892 0.21852 1.456759 1.115452
Tiempo que demora esta sección: 0 minutos
Generamos un gráfico de PPOO por cada conglomerado.
Código
ppoo_clus_pre_pam_om4_q<- df_filled[,c("run","glosa_pueblo_originario")] %>% dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om4","factor_inclusivo_real_hist_mas_autperc")], by="run", multiple="first") %>% dplyr::mutate(glosa_pueblo_originario_rec= dplyr::case_when(glosa_pueblo_originario=="NINGUNO"& factor_inclusivo_real_hist_mas_autperc!="00"~"DESCONOCIDO", T~glosa_pueblo_originario)) %>% janitor::tabyl(glosa_pueblo_originario_rec, clus_pam_om4) %>% janitor::adorn_percentages("row")reshape2::melt(ppoo_clus_pre_pam_om4_q, id.vars ="glosa_pueblo_originario_rec") %>% dplyr::mutate(glosa_pueblo_originario_rec= dplyr::recode(glosa_pueblo_originario_rec, "OTRO (ESPECIFICAR)"="OTRO(n=77)", "RAPA NUI (PASCUENSE)"="RAPA NUI(n=34)", "YAGÁN (YÁMANA)"="YAGÁN(n=2)","AYMARA"="AYMARA(n=13)","COLLA"="COLLA(n=6)","DIAGUITA"="DIAGUITA(n=3)","KAWÉSQAR"="KAWÉSQAR(n=4)","MAPUCHE"="MAPUCHE(n=255)","DESCONOCIDO"=".DESCONOCIDO(n=1.985)","NINGUNO"=".NINGUNO(n=9.156)")) %>%ggplot(aes(x = glosa_pueblo_originario_rec, y = value, fill = variable)) +geom_bar(stat ="identity", position ="fill") +scale_fill_manual(values =c("#D2B48C", "#E27A5B", "#708090", "#6B8E23")) +labs(title =NULL,x ="Grupo Étnico",y ="Proporción de Casos",fill ="Grupos") +# Cambia el título de la leyenda a "Grupos"theme_minimal() +theme(axis.text.y =element_text(size =12), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =12), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =14), # Tamaño del título del eje Xaxis.title.y =element_text(size =14), # Tamaño del título del eje Yplot.title =NULL, # Tamaño y estilo del título del gráficolegend.title =element_text(size =14, margin =margin(b =-.1)), # Tamaño del título de la leyendalegend.spacing.y =unit(1.5, "lines"),legend.box.spacing =unit(0.5, "lines"), # Controla el espacio entre la leyenda y el gráficolegend.margin =margin(5, 5, 5, 5), legend.key.height =unit(1, "cm"), legend.text =element_text(size =12) # Tamaño del texto de la leyenda ) +coord_flip() # Hacer el gráfico horizontalggsave("_figs/grafico_ancho_achatado_pam_om4_q.png", width =10, height =5, dpi=1000)ppoo_clus_pre_pam_om4_q_rapanui<- df_filled[,c("run","glosa_pueblo_originario")] %>% dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om4","factor_inclusivo_real_hist_mas_autperc")], by="run", multiple="first") %>% dplyr::mutate(glosa_pueblo_originario_rec= dplyr::case_when(glosa_pueblo_originario=="NINGUNO"& factor_inclusivo_real_hist_mas_autperc!="00"~"DESCONOCIDO", glosa_pueblo_originario=="NINGUNO"~"NINGUNO", glosa_pueblo_originario=="RAPA NUI (PASCUENSE)"~"RN", T~"RESTO")) %>% janitor::tabyl(glosa_pueblo_originario_rec, clus_pam_om4) %>% janitor::adorn_percentages("row")cat("origen?")df_filled[,c("run","glosa_pueblo_originario")] %>% dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run", "clus_pam_om4","factor_inclusivo_real_hist_mas_autperc", "codigo_region_rec_base")], by="run", multiple="first") %>% dplyr::mutate(glosa_pueblo_originario_rec= dplyr::case_when(glosa_pueblo_originario=="NINGUNO"& factor_inclusivo_real_hist_mas_autperc!="00"~"DESCONOCIDO", glosa_pueblo_originario=="NINGUNO"~"NINGUNO", glosa_pueblo_originario=="RAPA NUI (PASCUENSE)"~"RN", T~"RESTO")) %>% dplyr::filter(glosa_pueblo_originario_rec=="RN") |> janitor::tabyl(codigo_region_rec_base)
origen? codigo_region_rec_base n percent
RM 19 0.5588235
noRM 15 0.4411765
PPOO por cluster
Tiempo que demora esta sección: 0 minutos
1.1.1. Trayectorias
Vemos los gráficos de las trayectorias
Código
categories_pam_om4_q<-attr(States_Wide.seq_quarter_t_prim_adm_cens, "labels")new_labels <- categories_pam_om4_qnew_labels[which(categories_pam_om4_q =="Otras causas")] <-"Otras\ncausas"#new_labels[which(categories == "Consumo\nde sustancias")] <- "Consumo de\nsustancias"# Creamos un vector con las columnas llenando con NA si faltan valoressil_pam_om_clus4_q <-wcSilhouetteObs(as.dist(dist_quarter_om), pamRange_quarter_om$clustering$cluster4, measure="ASW")seq_plot_pam_om4_q <-ggseqiplot(States_Wide.seq_quarter_t_prim_adm_cens, group= ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4,facet_ncol=2, facet_nrow=2, sortv=sil_pam_om_clus4_q) +theme(legend.position ="none")+labs(x="Trimestres", y="# IDs de usuarios")+#guides(fill = guide_legend(nrow = 1))+theme(panel.spacing =unit(0.1, "lines"), # Reduce el espaciado entre los panelesaxis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-10)),#,margin = margin(l = -10)),strip.text =element_text(size =11, margin =margin(b =-15)),legend.text =element_text(size =15),legend.spacing.x =unit(0.1, 'cm'), # Alinea el título de la leyenda hacia la izquierdalegend.box.margin =margin(t =0, r =0, b =0, l =-50),legend.position ="bottom", legend.justification ="left",panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda )+guides(fill =guide_legend(nrow =1)) +scale_fill_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))+scale_color_manual(labels = new_labels, values=c("#E2725B", "#556B2F", "#D2B48C",#"#8B4513","#FFFFFF","#808080","#000000"))seq_plot_pam_om4_q ggsave(filename="_figs/clusters_pam_om4_q_mod.png", seq_plot_pam_om4_q, width =9.5, height =5.5, dpi=1000)
Trayectorias de hospitalización, orden de sujetos según el primer estado observado y su duración, representando a cada individuo como una línea en el gráfico (observaciones ordenadas de acuerdo a ASW)
Tiempo que demora esta sección: 0.3 minutos
Código
seq_plot2_pam_om4_q <-ggseqdplot(States_Wide.seq_quarter_t_prim_adm_cens, group= ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4,facet_ncol=2, facet_nrow=2) +theme(legend.position ="none")+# Colocar la leyenda abajolabs(x="Trimestres", y="Frecuencia relativa de estados")+theme(panel.spacing =unit(0.1, "lines"),axis.text.y =element_text(size =15), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =15), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =15), # Tamaño del título del eje Xaxis.title.y =element_text(size =15, margin =margin(r =-5)),strip.text =element_text(size =11),panel.spacing.y =unit(0.5, "lines"),strip.placement ="outside", # Para colocar las tiras fuera de los ejesstrip.background =element_blank() # Elimina el fondo para que parezca más espacioso#legend.key.size = unit(1.5, "lines"), # Aumenta el tamaño de los símbolos en la leyenda ) # Colocar la leyenda abajoseq_plot2_pam_om4_qggsave("_figs/clusterspam_om42_q_mod.png",seq_plot2_pam_om4_q, width =8.5, height =5.5, dpi=1000)table_data_pam_om4_q <-sprintf("%1.2f",pamRange_quarter_om$stats[3,])table_data_pam_om4_q <-as.data.frame(t(table_data_pam_om4_q))colnames(table_data_pam_om4_q)<-attr(pamRange_quarter_om$stats, "name")table_data_pam_om4_q %>% knitr::kable()
PBC
HG
HGSD
ASW
ASWw
CH
R2
CHsq
R2sq
HC
0.50
0.67
0.67
0.60
0.60
922.92
0.31
1269.60
0.39
0.18
Trayectorias de hospitalización, frecuencia relativa de estados en un gráfico de barras apiladas por trimestre.
Trayectorias de hospitalización, frecuencia relativa de estados en un gráfico de barras apiladas por trimestre.
Tiempo que demora esta sección: 0.1 minutos
De este modo, presenta el cambio agregado en la distribución de estados a lo largo del tiempo, sin considerar las secuencias individuales.
Código
invisible("Definimos las observaciones que tienen siluetas negativas")sil_neg_pam_om_clus4_q <-which(sil_pam_om_clus4_q<0)invisible("A qué conglomerados pertenecen?")table(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[sil_neg_pam_om_clus4_q, "clus_pam_om4"])ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$rn<-1:nrow(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens)
clus_pam_om4
6035, Un trimestre, TSM(4) 6025, Un trimestre, TUS(3) 5939, Un semestre TSM(1) 5989, Comorbilidad un trimestre(2)
2 0 133 0
Tiempo que demora esta sección: 0 minutos
1.1.2.Exploración transiciones
1.1.2.a Transiciones- RM y no RM
Tasas de transición no RM a RM y viceversa
Código
invisible("Tasas de transición no RM a RM y viceversa")trim_tasa_pam_om4_q_cens_cnt<-seqcount_t(States_Wide.seq_quarter_t_prim_adm_RM_cens, group=ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4) %>% dplyr::filter(count>0) %>% dplyr::mutate(trans =paste0(from,"_", to)) %>% dplyr::mutate(across(c("from","to"),~gsub("\\[->\\s*|\\s*->\\s*\\]|\\[|\\]", "", .))) trim_tasa_pam_om4_q_cens_rate<-seqtrate_t(States_Wide.seq_quarter_t_prim_adm_RM_cens, group=ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4) %>% dplyr::filter(rate>0) %>% dplyr::mutate(trans =paste0(from,"_", to)) %>% dplyr::mutate(across(c("from","to"),~gsub("\\[->\\s*|\\s*->\\s*\\]|\\[|\\]", "", .)))
Tiempo que demora esta sección: 0 minutos
Código
trim_tasa_pam_om4_q_cens_rate %>% dplyr::left_join(trim_tasa_pam_om4_q_cens_cnt, by=c("from"="from", "glosa_sexo"="glosa_sexo","to"="to")) %>% dplyr::rename("recuento"="count") %>% dplyr::filter(from %in%c("RM", "noRM")) %>%ggplot(aes(x = from, y = to, fill = rate, size=log(recuento+1))) +geom_tile() +coord_flip()+scale_fill_gradient(low ="white", high ="blue") +# Ajusta la escala de colores según tus preferenciaslabs(title ="Tasas de transición, Trimestre (s/censura)",x ="Desde",y ="Hacia",fill ="Rate") +theme_minimal() +facet_wrap(~glosa_sexo)+theme(axis.text.x =element_text(angle =45, hjust =1))+geom_text(aes(label =sprintf("%1.2f", rate), size =log(recuento+1)*.5), color ="black")invisible("Hay muy pocos casos que se entrecruzan entre noRM y RM (fuera de la diagnonal)")
Porcentajes de transición no-RM y RM por cada cluster
Tiempo que demora esta sección: 0.1 minutos
Hay muy pocos casos que se entrecruzan entre noRM y RM (fuera de la diagnonal)
trim_tasa2_pam_om4_q_cens_rate %>% dplyr::left_join(trim_tasa2_pam_om4_q_cens_cnt, by=c("from"="from", "glosa_sexo"="glosa_sexo","to"="to")) %>% dplyr::rename("recuento"="count") %>%#dplyr::filter(from %in% c("RM", "noRM")) %>% ggplot(aes(x = from, y = to, fill = rate, size=log(recuento+1))) +geom_tile() +coord_flip()+scale_fill_gradient(low ="white", high ="blue") +# Ajusta la escala de colores según tus preferenciaslabs(title ="Tasas de transición, Trimestre (s/censura)",x ="Desde",y ="Hacia",fill ="Rate") +theme_minimal() +facet_wrap(~glosa_sexo)+theme(axis.text.x =element_text(angle =45, hjust =1))+geom_text(aes(label =sprintf("%1.2f", rate), size =log(recuento+1)*.5), color ="black")
Porcentajes de transición, transiciones posteriores, por cada cluster
Tiempo que demora esta sección: 0.1 minutos
1.1.2.c Tiempo promedio por cluster
Código
seq_mean_t(States_Wide.seq_quarter_t_prim_adm_cens, group=ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4) %>% data.table::as.data.table(keep.rowname=T) %>% dplyr::mutate(rn=gsub("\\d", "", rn)) |>ggplot(aes(x=rn, fill= factor_inclusivo, y=Mean))+geom_bar(width =1, stat ="identity") +theme_minimal() +facet_wrap(~factor_inclusivo)+labs(title =NULL,x =NULL,y =NULL) +scale_fill_manual(values =c("#70809090", "#6B8E2380", "#E27A5B","#D2B48C")) +coord_flip()+theme(#axis.text.x = element_blank(),#axis.text.y = element_blank(),panel.grid =element_blank()) +# scale_fill_brewer(palette = "Pastel1", labels=c("Sin\nautoidentificación\nni reconocimiento", "Autoidentificación\nsin reconocimiento", "Ambas")) +geom_text(aes(label =round(Mean,1)), position =position_stack(vjust =0.5), size =3.5, # Ajusta el tamaño de la fuente aquícolor ="black", # Color del textofamily ="sans", # Puedes cambiar la fuente si lo deseasbackground =element_rect(fill ="white", color =NA)) +# Fondo blancotheme(legend.title =element_blank())invisible("No me aporta mucho")
Tiempo promedio en cada estado por estatus PPOO (Trimestral c/censura)
Tiempo que demora esta sección: 0.1 minutos
Observamos que aquellos en el conglomerado que se encuentra un trimestre, la duración de los ingresos relacionados con trastornos de salud mental son en promedio más largos.
1.1.3. Comparación variables
1.1.3.a. Comparación covariables- PPOO
Código
ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>% dplyr::count(clus_pam_om4, factor_inclusivo_real_hist_mas_autperc) %>% dplyr::group_by(clus_pam_om4) %>% dplyr::mutate(n_prop =paste0(n, " (",scales::percent(n /sum(n), accuracy=.1),")")) %>% dplyr::select(-n) %>% tidyr::pivot_wider(names_from = factor_inclusivo_real_hist_mas_autperc, values_from = n_prop, values_fill ="0") %>% knitr::kable("markdown", col.names=c("Conglomerados","No se identifica/no pertenece", "No se identifica/hay reconocimiento", "Se identifica/hay reconocimeinto"), caption="Porcentajes por fila, conglomerado vs. Pertenencia/identificación + Reconocimento CONADI PPOO")invisible("6025 tiwnw un poxo mas PPOO, lo mismo con 5710")
Porcentajes por fila, conglomerado vs. Pertenencia/identificación + Reconocimento CONADI PPOO
Conglomerados
No se identifica/no pertenece
No se identifica/hay reconocimiento
Se identifica/hay reconocimeinto
6035, Un trimestre, TSM(4)
3857 (80.6%)
536 (11.2%)
395 (8.2%)
6025, Un trimestre, TUS(3)
537 (78.5%)
79 (11.5%)
68 (9.9%)
5939, Un semestre TSM(1)
291 (82.2%)
39 (11.0%)
24 (6.8%)
5989, Comorbilidad un trimestre(2)
178 (84.0%)
18 (8.5%)
16 (7.5%)
Tiempo que demora esta sección: 0 minutos
Vemos las categorías de clasificación de PPOO según autopercepción (en MINSAL y en RSH) y reconocimiento CONADI.
ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>% dplyr::count(clus_pam_om4, factor_inclusivo_real_hist_mas_autperc_bin) %>% dplyr::group_by(clus_pam_om4) %>% dplyr::mutate(n_prop =paste0(n, " (",scales::percent(n /sum(n), accuracy=.1),")")) %>% dplyr::select(-n) %>% tidyr::pivot_wider(names_from = factor_inclusivo_real_hist_mas_autperc_bin, values_from = n_prop, values_fill ="0") %>% knitr::kable("markdown", col.names=c("Conglomerados","No se identifica/no pertenece", "Hay reconocimeinto"), caption="Porcentajes por fila, conglomerado vs. Pertenencia/identificación + Reconocimento CONADI PPOO")invisible("6025 tiwnw un poxo mas PPOO, lo mismo con 5710")
Porcentajes por fila, conglomerado vs. Pertenencia/identificación + Reconocimento CONADI PPOO
Conglomerados
No se identifica/no pertenece
Hay reconocimeinto
6035, Un trimestre, TSM(4)
3857 (80.6%)
931 (19.4%)
6025, Un trimestre, TUS(3)
537 (78.5%)
147 (21.5%)
5939, Un semestre TSM(1)
291 (82.2%)
63 (17.8%)
5989, Comorbilidad un trimestre(2)
178 (84.0%)
34 (16.0%)
Tiempo que demora esta sección: 0 minutos
1.1.3.b. Comparación covariables- Mortalidad
Código
# invisible("No hay nada, el tiempo promedio de censura es similar")ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>% dplyr::mutate(death_time_rec=ifelse(death_time==20,0,1)) %>% janitor::tabyl(clus_pam_om4,death_time_rec) %>% dplyr::mutate(`1`=paste0(`1`," (", scales::percent(`1`/(`0`+`1`), accuracy=.1),")")) %>% dplyr::left_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>% dplyr::group_by(clus_pam_om4) %>% dplyr::summarise(mean=sprintf("%1.1f",mean(cens_time))), by="clus_pam_om4") %>% dplyr::select(-`0`) %>% knitr::kable("markdown", col.names=c("Conglomerado","Mortalidad observada", "Promedio"), caption="Post-hoc, conglomerado vs. Mortalidad y tiempo a censura")
Post-hoc, conglomerado vs. Mortalidad y tiempo a censura
Conglomerado
Mortalidad observada
Promedio
6035, Un trimestre, TSM(4)
45 (0.9%)
17.9
6025, Un trimestre, TUS(3)
13 (1.9%)
18.2
5939, Un semestre TSM(1)
8 (2.3%)
18.0
5989, Comorbilidad un trimestre(2)
4 (1.9%)
18.1
Tiempo que demora esta sección: 0 minutos
Código
ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>% dplyr::mutate(death_time_rec=ifelse(death_time==20,0,1)) %>% janitor::tabyl(death_time_rec,clus_pam_om4) %>% janitor::chisq.test(correct=T)#X-squared = 10.014, df = 3, p-value = 0.01845fisher.test(with(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens%>% dplyr::mutate(death_time_rec=ifelse(death_time==20,0,1)), table(death_time_rec , clus_pam_om4), simulate.p.value=T, B=1e5))message("Descartando valores negativos en sil width")chisq_cramerv(with(subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, !rn %in% sil_neg_pam_om_clus4_q)%>% dplyr::mutate(death_time_rec=ifelse(death_time==20,0,1)), table(death_time_rec , clus_pam_om4)))fisher.test(with(subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, !rn %in% sil_neg_pam_om_clus4_q)%>% dplyr::mutate(death_time_rec=ifelse(death_time==20,0,1)), table(death_time_rec , clus_pam_om4), simulate.p.value=T, B=1e5))invisible("no se basa en la distribución chi-cuadrado. Fisher se basa en permutaciones exactas, por lo que no se calculan df.")#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_##_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_#_tab_cl_mortalidad_pam_om4_q<- ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>% dplyr::mutate(death_time_rec=ifelse(death_time==20,0,1)) |> janitor::tabyl(death_time_rec,clus_pam_om4) |>as.matrix(ncol=2)labels_pam_om4_q <-c("6035, Un trimestre, TSM(4)","6025, Un trimestre, TUS(3)","5939, Un semestre TSM(1)","5989, Comorbilidad un trimestre(2)")# Realizar el análisis y crear la tabla directamentepairwise.prop.test(t(tab_cl_mortalidad_pam_om4_q[,2:5]), p.adjust.method ="holm")$p.value |>as.table() |>as.data.frame() |>rename(Grupo_1 = Var1, Grupo_2 = Var2, p_value = Freq) |>filter(!is.na(p_value)) |>mutate(Grupo_1 = labels_pam_om4_q[as.numeric(Grupo_1)],Grupo_2 = labels_pam_om4_q[as.numeric(Grupo_2)],p_value =ifelse(p_value <.001, "<.001", sprintf("%1.3f",p_value)) ) |>kable(col.names =c("Grupo 1", "Grupo 2", "Valor p ajustado"),align ="l",caption="Corrección parcial por comparaciones múltiples (Holm–Bonferroni)" )
Pearson's Chi-squared test
data: .
X-squared = 10.014, df = 3, p-value = 0.01845
Fisher's Exact Test for Count Data
data: with(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>% dplyr::mutate(death_time_rec = ifelse(death_time == 20, 0, 1)), table(death_time_rec, clus_pam_om4), simulate.p.value = T, B = 1e+05)
p-value = 0.01251
alternative hypothesis: two.sided
$chisq_statistic
[1] "8.76"
$chisq_df
df
3
$chisq_p_value
[1] "0.0326"
$cramers_v
[1] "0.04"
Fisher's Exact Test for Count Data
data: with(subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, !rn %in% sil_neg_pam_om_clus4_q) %>% dplyr::mutate(death_time_rec = ifelse(death_time == 20, 0, 1)), table(death_time_rec, clus_pam_om4), simulate.p.value = T, B = 1e+05)
p-value = 0.02068
alternative hypothesis: two.sided
Corrección parcial por comparaciones múltiples (Holm–Bonferroni)
Grupo 1
Grupo 2
Valor p ajustado
6035, Un trimestre, TSM(4)
6035, Un trimestre, TSM(4)
0.214
6025, Un trimestre, TUS(3)
6035, Un trimestre, TSM(4)
0.214
5939, Un semestre TSM(1)
6035, Un trimestre, TSM(4)
1.000
6025, Un trimestre, TUS(3)
6025, Un trimestre, TUS(3)
1.000
5939, Un semestre TSM(1)
6025, Un trimestre, TUS(3)
1.000
5939, Un semestre TSM(1)
5939, Un semestre TSM(1)
1.000
Tiempo que demora esta sección: 0 minutos
Código
# Cargar las librerías necesariaslibrary(survival)
Adjuntando el paquete: ‘survival’
The following object is masked from ‘package:survminer’:
myeloma
Código
library(ggplot2)# Crear la variable de supervivenciaing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$surv_obj_4c <-Surv(time = ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$death_time,event =ifelse(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$death_time==20,0,1))# Realizar el análisis de Log-Rank (survdiff)surv_diff_4c <-survdiff(surv_obj_4c ~ clus_pam_om4,data = ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens)# Mostrar los resultados del testprint(surv_diff_4c)
Call:
survdiff(formula = surv_obj_4c ~ clus_pam_om4, data = ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens)
N Observed Expected (O-E)^2/E (O-E)^2/V
clus_pam_om4=6035, Un trimestre, TSM(4) 4788 45 55.54 2.000 9.68
clus_pam_om4=6025, Un trimestre, TUS(3) 684 13 7.91 3.273 3.69
clus_pam_om4=5939, Un semestre TSM(1) 354 8 4.10 3.702 3.93
clus_pam_om4=5989, Comorbilidad un trimestre(2) 212 4 2.45 0.988 1.02
Chisq= 10 on 3 degrees of freedom, p= 0.02
Pairwise comparisons using Log-Rank test
data: ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens and clus_pam_om4
6035, Un trimestre, TSM(4) 6025, Un trimestre, TUS(3) 5939, Un semestre TSM(1)
6025, Un trimestre, TUS(3) 0.11 - -
5939, Un semestre TSM(1) 0.11 1.00 -
5989, Comorbilidad un trimestre(2) 0.68 1.00 1.00
P value adjustment method: holm
Código
message("Sin siluetas negativas")
Sin siluetas negativas
Código
surv_diff_4c_neg_sil <-survdiff(Surv(time = death_time,event =ifelse(death_time==20,0,1)) ~ clus_pam_om4,data =subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, !rn %in% sil_neg_pam_om_clus4_q))# Mostrar los resultados del testprint(surv_diff_4c_neg_sil)
Call:
survdiff(formula = Surv(time = death_time, event = ifelse(death_time ==
20, 0, 1)) ~ clus_pam_om4, data = subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens,
!rn %in% sil_neg_pam_om_clus4_q))
N Observed Expected (O-E)^2/E (O-E)^2/V
clus_pam_om4=6035, Un trimestre, TSM(4) 4786 45 54.35 1.61 8.53
clus_pam_om4=6025, Un trimestre, TUS(3) 684 13 7.75 3.56 4.03
clus_pam_om4=5939, Un semestre TSM(1) 221 5 2.51 2.48 2.58
clus_pam_om4=5989, Comorbilidad un trimestre(2) 212 4 2.39 1.08 1.12
Chisq= 8.7 on 3 degrees of freedom, p= 0.03
Pairwise comparisons using Log-Rank test
data: subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, and clus_pam_om4 !rn %in% sil_neg_pam_om_clus4_q) and clus_pam_om4
6035, Un trimestre, TSM(4) 6025, Un trimestre, TUS(3) 5939, Un semestre TSM(1)
6025, Un trimestre, TUS(3) 0.13 - -
5939, Un semestre TSM(1) 0.27 1.00 -
5989, Comorbilidad un trimestre(2) 0.68 1.00 1.00
P value adjustment method: holm
Código
# Ajustar el modelo de Kaplan-Meierkm_fit_4c <-survfit(surv_obj_4c ~ clus_pam_om4,data = ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens)# Extraer los datos del modelo Kaplan-Meier para usar con ggplotkm_data_4c <-data.frame(time = km_fit_4c$time,surv = km_fit_4c$surv,upper = km_fit_4c$upper,lower = km_fit_4c$lower,strata =rep( c("6035,\nUn trimestre,\nTSM(4)", "6025,\nUn trimestre,\nTUS(3)", "5939,\nUn semestre\nTSM(1)", "5989,\nComorbilidad un\ntrimestre(2)"), km_fit_4c$strata))biostat3::survRate(Surv(time = death_time,event =ifelse(death_time==20,0,1)) ~ clus_pam_om4, data= ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens) %>% dplyr::mutate(across(c("rate", "lower", "upper"),~sprintf("%1.2f",.*10000)))
clus_pam_om4 tstop event rate lower upper
clus_pam_om4=6035, Un trimestre, TSM(4) 6035, Un trimestre, TSM(4) 95127.012 45 4.73 3.45 6.33
clus_pam_om4=6025, Un trimestre, TUS(3) 6025, Un trimestre, TUS(3) 13518.762 13 9.62 5.12 16.44
clus_pam_om4=5939, Un semestre TSM(1) 5939, Un semestre TSM(1) 7002.082 8 11.43 4.93 22.51
clus_pam_om4=5989, Comorbilidad un trimestre(2) 5989, Comorbilidad un trimestre(2) 4181.071 4 9.57 2.61 24.50
Código
# Crear el gráfico de Kaplan-Meier con ggplot2ggplot(km_data_4c, aes(x = time, y = surv, color = strata)) +geom_step(size =1.2) +# Curvas de supervivencia#geom_ribbon(aes(ymin = lower, ymax = upper, fill = strata), alpha = 0.2, color = NA) + # Intervalos de confianzalabs(title ="Curvas de Kaplan-Meier",x ="Tiempo (meses)",y ="Probabilidad de Supervivencia",color ="Grupo",fill ="Grupo" ) +theme_minimal() +theme(legend.position ="bottom")
Al parecer, el conglomerado 6025, Un trimestre, TUS(3) tendría una proporción significativamente mayor que el resto de pacientes fuera de la RM, aunque con una fuerza de asociación débil. En contraste, el conglomerado 5989, Comorbilidad un trimestre(2) tiene una proporción mayor de personas en la RM a diferencia de 6025, Un trimestre, TUS(3) y 6035, Un trimestre, TSM(4).
Código
pairwise_chisq_gof_test(table(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$clus_pam_om4,ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens$codigo_region_rec_base), p.adjust.method="holm")|> knitr::kable("markdown", caption="Dependencia categórica sol. 4 conglomerados, por pares de categorías en RM")
Dependencia categórica sol. 4 conglomerados, por pares de categorías en RM
De la tabla anterior, se observa que participantes pertenecientes al conglomerado “6025, Un trimestre, TUS(3)” presenta una menor proporción de residentes de la Región Metropolitana. Por otra parte, participantes pertenecientes al conglomerado ” 5989, Comorbilidad un trimestre(2)” pertenecen en mayor proporción a la RM.
1.1.3.e. Comparación covariables- Región
Código
tab_cluster_region_pam_om4_q<-ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens %>% dplyr::inner_join(data_long_establecimiento_2024_std[,c("ESTAB_HOMO", "codigo_region", "nivel_de_atencion", "nivel_de_complejidad")], by =c("estab_homo_base"="ESTAB_HOMO"), multiple ="first") %>% janitor::tabyl(codigo_region, clus_pam_om4) %>% janitor::adorn_percentages("col") %>% janitor::adorn_rounding(digits =2)#colnames(tab_cluster_region_pam_om4_q)<- c("reg", "c1", "c4", "c3", "c5", "c6", "c7", "c8", "c9", "c2")cod_reg_homo_pam_om4_q<-data.frame(codigo_region =1:16,nombre_region =c("Región de Tarapacá","Región de Antofagasta","Región de Atacama","Región de Coquimbo","Región de Valparaíso","Región del Libertador General Bernardo O'Higgins","Región del Maule","Región del Biobío","Región de La Araucanía","Región de Los Lagos","Región de Aysén del General Carlos Ibáñez del Campo","Región de Magallanes y de la Antártica Chilena","Región Metropolitana de Santiago","Región de Los Ríos","Región de Arica y Parinacota","Región de Ñuble" ),stringsAsFactors =FALSE)dplyr::mutate(tab_cluster_region_pam_om4_q, promedio_fila =rowMeans(across(2:length(colnames(tab_cluster_region_pam_om4_q))))) %>% dplyr::arrange(desc(promedio_fila)) %>% dplyr::left_join(cod_reg_homo_pam_om4_q, by="codigo_region") %>% dplyr::select(codigo_region, nombre_region, everything()) %>% dplyr::select(-promedio_fila) %>% dplyr::mutate_at(3:(length(colnames(tab_cluster_region_pam_om4_q))+1),~scales::percent(.)) %>% knitr::kable(caption="Porcentaje por región")
Porcentaje por región
codigo_region
nombre_region
6035, Un trimestre, TSM(4)
6025, Un trimestre, TUS(3)
5939, Un semestre TSM(1)
5989, Comorbilidad un trimestre(2)
13
Región Metropolitana de Santiago
45%
35%
47%
56%
10
Región de Los Lagos
6%
20%
5%
10%
5
Región de Valparaíso
8%
13%
11%
6%
8
Región del Biobío
10%
8%
9%
9%
9
Región de La Araucanía
5%
4%
5%
4%
6
Región del Libertador General Bernardo O’Higgins
4%
3%
5%
2%
7
Región del Maule
4%
3%
3%
4%
2
Región de Antofagasta
2%
1%
4%
1%
11
Región de Aysén del General Carlos Ibáñez del Campo
pairwise_chisq_gof_test(dplyr::filter(tab_clus_macrozona_pam_om4_q,macrozona!="RM")[-1], p.adjust.method="holm")|> knitr::kable("markdown", caption="Dependencia categórica sol. 4 conglomerados, por pares de categorías en Macrozona (corrección Holm-Bonferroni)")#Groups sharing a letter are not significantlt different (alpha = 0.05)
Dependencia categórica sol. 4 conglomerados, por pares de categorías en Macrozona (corrección Holm-Bonferroni)
n
group1
group2
statistic
p
df
p.adj
p.adj.signif
3082
6035, Un trimestre, TSM(4)
6025, Un trimestre, TUS(3)
70.575017
0.00e+00
4
0.00e+00
****
2829
6035, Un trimestre, TSM(4)
5939, Un semestre TSM(1)
3.767090
4.38e-01
4
4.38e-01
ns
2733
6035, Un trimestre, TSM(4)
5989, Comorbilidad un trimestre(2)
8.456608
7.62e-02
4
2.29e-01
ns
631
6025, Un trimestre, TUS(3)
5939, Un semestre TSM(1)
30.008475
4.90e-06
4
2.44e-05
****
535
6025, Un trimestre, TUS(3)
5989, Comorbilidad un trimestre(2)
6.949579
1.39e-01
4
2.78e-01
ns
282
5939, Un semestre TSM(1)
5989, Comorbilidad un trimestre(2)
10.498374
3.28e-02
4
1.31e-01
ns
Tiempo que demora esta sección: 0 minutos
Hay una asociación entre pacientes en 6025, Un trimestre, TUS(3) y el resto de conglomerados. Particularmente, pacientes de este conglomerado pertenecen en mayor proporción a la zona sur (La Araucanía, Los Ríos y Los Lagos) y centro (Región de Coquimbo y de Valparaíso), aunque menos a la macrozona norte (Arica y Parinacota, Tarapacá, Antofagasta y Atacama) y la RM. Hay una asociación entre Pacientes en 6035, Un trimestre, TSM(4), en que éstos parecen pertenecer en mayor proporción a la Macrozona Norte, a excepción de comparado con pacientes en 5939, Un semestre TSM(1). Por otra parte, pacientes en el conglomerado 5989, Comorbilidad un trimestre(2) se encuentran en mayor proporción en la región Metropolitana
pairwise_chisq_gof_test(tab_clus_sexo_pam_om4_q[-1], p.adjust.method="holm")|> dplyr::mutate(p.adj= dplyr::case_when(p.adj<0.001~"<0.001",T~sprintf("%1.3f",p.adj)))|> knitr::kable("markdown", caption="Dependencia categórica sol. 4 conglomerados, por pares de categorías en Sexo (corrección Holm-Bonferroni)")
Dependencia categórica sol. 4 conglomerados, por pares de categorías en Sexo (corrección Holm-Bonferroni)
n
group1
group2
statistic
p
df
p.adj
p.adj.signif
5472
6035, Un trimestre, TSM(4)
6025, Un trimestre, TUS(3)
195.8119436
0.000
1
<0.001
****
5142
6035, Un trimestre, TSM(4)
5939, Un semestre TSM(1)
0.0356785
0.850
1
1.000
ns
5000
6035, Un trimestre, TSM(4)
5989, Comorbilidad un trimestre(2)
60.4441006
0.000
1
<0.001
****
1038
6025, Un trimestre, TUS(3)
5939, Un semestre TSM(1)
81.8920873
0.000
1
<0.001
****
896
6025, Un trimestre, TUS(3)
5989, Comorbilidad un trimestre(2)
0.0766505
0.782
1
1.000
ns
566
5939, Un semestre TSM(1)
5989, Comorbilidad un trimestre(2)
40.2303954
0.000
1
<0.001
****
Tiempo que demora esta sección: 0 minutos
Pacientes en 6025, Un trimestre, TUS(3) y 5989, Comorbilidad un trimestre(2) están compuesto por menos mujeres que 6035, Un trimestre, TSM(4) y 5939, Un semestre TSM(1).
Edad promedio primer ingreso con intervalo de confianza por conglomerado
Tiempo que demora esta sección: 0.1 minutos
Código
invisible("Prueba de Levene par igualdad de varianzas")with(dt_ing_calendar_quarter_t_desde_primera_adm_dedup %>% dplyr::filter(quarter ==0) %>% dplyr::inner_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run","clus_pam_om4")], by="run"), car::leveneTest(min_edad_anos, clus_pam_om4))#NSanova_clus_pam_om4_q <-oneway.test(min_edad_anos ~ clus_pam_om4, data = dt_ing_calendar_quarter_t_desde_primera_adm_dedup %>% dplyr::filter(quarter ==0) %>% dplyr::inner_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run","clus_pam_om4")], by="run"),var.equal = T)with(dt_ing_calendar_quarter_t_desde_primera_adm_dedup %>% dplyr::filter(quarter ==0) %>% dplyr::inner_join(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens[,c("run","clus_pam_om4")], by="run"), oneway_anova_effect_size (min_edad_anos, clus_pam_om4))# Ver los resultados del ANOVAprint(anova_clus_pam_om4_q)#F = 51.5, num df = 3, denom df = 6034, p-value < 2.2e-16#$eta_squared#[1] 0.0249659message("Descartando valores negativos en sil width")with(dt_ing_calendar_quarter_t_desde_primera_adm_dedup %>% dplyr::filter(quarter ==0) %>% dplyr::inner_join(subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, !rn %in% sil_neg_pam_om_clus4_q)[,c("run","clus_pam_om4")], by="run"), oneway.test(min_edad_anos ~ clus_pam_om4,var.equal = T) )
Levene's Test for Homogeneity of Variance (center = median)
Df F value Pr(>F)
group 3 1.2631 0.2852
6034
$anova_summary
Df Sum Sq Mean Sq F value Pr(>F)
group 3 2853 951.1 51.5 <2e-16 ***
Residuals 6034 111440 18.5
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
$eta_squared
[1] 0.0249659
One-way analysis of means
data: min_edad_anos and clus_pam_om4
F = 51.5, num df = 3, denom df = 6034, p-value < 2.2e-16
One-way analysis of means
data: min_edad_anos and clus_pam_om4
F = 49.735, num df = 3, denom df = 5899, p-value < 2.2e-16
Post-hoc, conglomerado vs. Promedio días de tratamiento
Conglomerado1
Conglomerado2
Estimación
6035, Un trimestre, TSM(4)
6025, Un trimestre, TUS(3)
1.96 [1.50, 2.41], p= <0.001
6035, Un trimestre, TSM(4)
5939, Un semestre TSM(1)
-0.82 [-1.41, -0.22], p= 0.002
6035, Un trimestre, TSM(4)
5989, Comorbilidad un trimestre(2)
1.14 [0.37, 1.90], p= <0.001
6025, Un trimestre, TUS(3)
5939, Un semestre TSM(1)
-2.77 [-3.48, -2.06], p= <0.001
6025, Un trimestre, TUS(3)
5989, Comorbilidad un trimestre(2)
-0.82 [-1.68, 0.04], p= 0.066
5939, Un semestre TSM(1)
5989, Comorbilidad un trimestre(2)
1.95 [1.02, 2.89], p= <0.001
Tiempo que demora esta sección: 0 minutos
Si bien la asociación entre edad mínima de ingreso es signifciatva, tiene un tamaño del efecto débil. Particularmente, pacientes clasificados como 6025, Un trimestre, TUS(3) y 5989, Comorbilidad un trimestre(2) ingresan de manera más tardía a tratamiento, a diferencia de pacientes clasificados en 5939, Un semestre TSM(1), quienes son los que ingresan más temprano, seguidos por una leve aunque significativa diferencia con 6035, Un trimestre, TSM(4).
Pearson's Chi-squared test
data: .
X-squared = 25.206, df = 12, p-value = 0.01388
Fisher's Exact Test for Count Data with simulated p-value (based on 1e+05 replicates)
data: .
p-value = 0.01067
alternative hypothesis: two.sided
Porcentajes por columna, conglomerado vs. Previsión
$chisq_statistic
[1] "25.21"
$chisq_df
df
12
$chisq_p_value
[1] "0.0139"
$cramers_v
[1] "0.04"
Fisher's Exact Test for Count Data with simulated p-value (based on 1e+05 replicates)
data: with(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, table(prev_benef_rec_post, clus_pam_om4))
p-value = 0.0107
alternative hypothesis: two.sided
$chisq_statistic
[1] "25.18"
$chisq_df
df
12
$chisq_p_value
[1] "0.0140"
$cramers_v
[1] "0.04"
Fisher's Exact Test for Count Data with simulated p-value (based on 1e+05 replicates)
data: with(subset(ing_dt_ing_calendar_quarter_t_desde_primera_adm_dedup_wide2_cens, !rn %in% sil_neg_pam_om_clus4_q), table(prev_benef_rec_post, clus_pam_om4))
p-value = 0.00993
alternative hypothesis: two.sided
Para post-hoc convendría consultar https://adaptivedesignstrial.shinyapps.io/posthoc/
Comparación post-hoc, conglomerado-previsión
Dimension
6035, Un trimestre, TSM(4)
6025, Un trimestre, TUS(3)
5939, Un semestre TSM(1)
5989, Comorbilidad un trimestre(2)
FFAA
3.25 (p=0.023)
-2.77 (p=0.112)
-0.42 (p=1.000)
-1.85 (p=1.000)
FONASA A
-1.69 (p=1.000)
0.78 (p=1.000)
1.73 (p=1.000)
0.18 (p=1.000)
FONASA BC
-0.52 (p=1.000)
-0.21 (p=1.000)
-0.24 (p=1.000)
1.80 (p=1.000)
FONASA D
-1.44 (p=1.000)
2.17 (p=0.601)
-0.72 (p=1.000)
0.36 (p=1.000)
ISAPRE
1.94 (p=1.000)
-0.97 (p=1.000)
-0.78 (p=1.000)
-1.58 (p=1.000)
Dependencia categórica sol. 4 conglomerados, por pares de categorías en Previsión (corrección Holm-Bonferroni)
n
group1
group2
statistic
p
df
p.adj
p.adj.signif
5472
6035, Un trimestre, TSM(4)
6025, Un trimestre, TUS(3)
13.913777
0.00758
4
0.0455
*
5142
6035, Un trimestre, TSM(4)
5939, Un semestre TSM(1)
3.934674
0.41500
4
1.0000
ns
5000
6035, Un trimestre, TSM(4)
5989, Comorbilidad un trimestre(2)
8.474783
0.07570
4
0.3780
ns
1038
6025, Un trimestre, TUS(3)
5939, Un semestre TSM(1)
5.224236
0.26500
4
1.0000
ns
896
6025, Un trimestre, TUS(3)
5989, Comorbilidad un trimestre(2)
3.040659
0.55100
4
1.0000
ns
566
5939, Un semestre TSM(1)
5989, Comorbilidad un trimestre(2)
4.970130
0.29000
4
1.0000
ns
Tiempo que demora esta sección: 0 minutos
Los pacientes en el conglomerado 6035, Un trimestre, TSM(4) muestra un porcentaje mayor de personas con previsión de las FFAA, particularmente, a diferencia de pacientes en el conglomerado 6025, Un trimestre, TUS(3).
Dependencia categórica sol. 4 conglomerados, por pares de categorías en Niv. Complejidad (corrección Holm-Bonferroni)
n
group1
group2
statistic
p
df
p.adj
p.adj.signif
5472
6035, Un trimestre, TSM(4)
6025, Un trimestre, TUS(3)
35.60142
3.00e-07
4
1.70e-06
****
5142
6035, Un trimestre, TSM(4)
5939, Un semestre TSM(1)
21.68587
2.31e-04
4
6.93e-04
***
5000
6035, Un trimestre, TSM(4)
5989, Comorbilidad un trimestre(2)
17.08293
1.86e-03
4
3.72e-03
**
1038
6025, Un trimestre, TUS(3)
5939, Un semestre TSM(1)
59.77247
0.00e+00
4
0.00e+00
****
896
6025, Un trimestre, TUS(3)
5989, Comorbilidad un trimestre(2)
14.72857
5.30e-03
4
5.30e-03
**
566
5939, Un semestre TSM(1)
5989, Comorbilidad un trimestre(2)
34.06543
7.00e-07
4
2.90e-06
****
Tiempo que demora esta sección: 0 minutos
Se constata una asociación significativa aunque débil entre nivel de complejidad de los establecimientos y el conglomerado en los que son clasificados los pacientes. Particularmente, quienes son clasificados en 5939, Un semestre TSM(1) se encuentran en menor proporción en centros de alta complejidad, mientras que pacientes en 5989, Comorbilidad un trimestre(2) se encuentran en mayor proporción en centros de alta complejidad. Por otra parte, pacientes clasificados en 6025, Un trimestre, TUS(3) se encuentran en menor proporción en centros de Mediana complejidad y en mayor proporción en centros de alta complejidad.
1.1.4. Compilación comparación covariables
Código
# Definir los datos correctamentedata_pam_om4_q <-cbind.data.frame(Grupo=c("6035,\nUn trimestre,\nTSM(4)", "6025,\nUn trimestre,\nTUS(3)", "5939,\nUn semestre\nTSM(1)", "5989,\nComorbilidad un\ntrimestre(2)"), PPOO_bin =c(NA, NA, NA, NA), PPOO_sinautoid =c(NA, NA, NA, NA), PPOO_conautoid =c(NA, NA, NA, NA), Mortalidad =c(NA, NA, NA, NA), RM =c(NA, "-", NA, "+"), `Macrozona-Austral`=c(NA, NA, NA, NA), `Macrozona-Centro`=c(NA, NA, NA, NA), `Macrozona-Centro Sur`=c(NA, NA, NA, NA), `Macrozona-Norte`=c("+","-", NA, NA), `Macrozona-Sur`=c(NA, "+", NA, NA), Sexo_mujeres =c(NA, "-", NA, "-"), `Edad ingreso`=c("-", "+", "-", "+"), `Previsión-FFAA`=c("+", "-", NA, NA), `Previsión-FONASA A`=c(NA, NA, NA, NA), `Previsión-FONASA BC`=c(NA, NA, NA, NA), `Previsión-FONASA D`=c(NA, NA, NA, NA), `Previsión-ISAPRE`=c(NA, NA, NA, NA), `NivComp-Baja`=c(NA, NA, NA, NA), `NivComp-Media`=c(NA, "-", NA, NA), `NivComp-Alta`=c(NA, "+", "-", "+"))## Asegurar que los nombres de las columnas sean válidos y no haya espacios en blanco# Derretir el dataframe para que sea adecuado para ggplot2data_melt_pam_om4_q <- reshape2::melt(data_pam_om4_q, id.vars ='Grupo', variable.name ='Variable', value.name ='Asociación')# Reemplazar los NA por un valor vacíodata_melt_pam_om4_q$Asociación[is.na(data_melt_pam_om4_q$Asociación)] <-"\n"# Crear el gráfico con ggplotplot_data_melt_pam_om4_q<-data_melt_pam_om4_q %>% dplyr::mutate(Variable =gsub("_", " ", Variable)) %>%ggplot(aes(x = Variable, y = Grupo, fill = Asociación)) +geom_tile(color ="white", size =0.8) +scale_fill_manual(values =c("+"="#556B2F", "-"="#E2725B", "\n"="white")) +labs(title =NULL, x ="Variables", y ="Conglomerado") +theme_minimal() +theme(#axis.text.x = element_text(angle = 45, hjust = 1),panel.grid =element_blank())+theme(axis.text.y =element_text(size =17*.7, face ="bold"),#,margin = margin(l = 7)), # Tamaño de las etiquetas de los grupos étnicosaxis.text.x =element_text(size =17*.7, face ="bold"), # Tamaño de las etiquetas del eje Xaxis.title.x =element_text(size =16*.7, face ="bold"),#,margin = margin(t = -15)), # Tamaño del título del eje Xaxis.title.y =element_text(size =16*.7, face ="bold"), # Tamaño del título del eje Yplot.title =NULL, # Tamaño y estilo del título del gráficolegend.title =element_text(size =17*.7, face ="bold"), # Tamaño del título de la leyendalegend.spacing.y =unit(1.2, "lines"),legend.box.spacing =unit(0.5, "lines"), # Controla el espacio entre la leyenda y el gráficolegend.margin =margin(5, 5, 5, 5), legend.key.height =unit(1, "cm"), legend.text =element_text(size =15*.7, face ="bold") # Tamaño del texto de la leyenda ) +coord_flip()plot_data_melt_pam_om4_qggsave(filename="_figs/asociaciones_pam_om4_q.png", plot_data_melt_pam_om4_q, width=11*.8, height=6*.8, dpi=600)
